//---------------------------------------------------------------------------
#include <stdio.h>
#include <string.h>
#include <limits.h>
#include <stdlib.h>
#ifdef BC
  #include <algorithm.h> // max min
#endif
//---------------------------------------------------------------------------
#include "main.h"
#include "conf.h"
#include "cover.h"
#include "boom.h"
#include "tboolfnct.h"
#include "tpfnct.h"
#include "decomp.h"
//---------------------------------------------------------------------------
void FindCover( TBoolFnct * fn, TGauge gauge )
{
  TCovs * c;
  TPfnct * ff;
  TPfnct * pff;
  int i, j, k, l;

  gauge.NewGauge( "Determining cover", fn->terms, 1 );
  ff = fn->fnct;
  l = 0;
  while ( ff != 0 )
  {
    ff->cov = 0;
    delete[] ff->covfield;
    ff->covlen = 0;
    ff->covfield = 0;
    for ( j = 0; j < fn->outvars; j++ )
      if ( ff->ov[j] == '1' )
      {
        pff = fn->parent->onset[j];
        k = 0;
        while ( pff != 0 )
        {
          for ( i = 0; i < fn->parent->inpvars; i++ )
            if ( ff->iv[i] != pff->iv[i] && ff->iv[i] != '-' ) break;
          if ( i == fn->parent->inpvars )
          {
            c = new TCovs;
            c->var = j;
            c->index = k;
            c->fn = pff;
            c->next = ff->cov;
            ff->cov = c;
            ff->covlen++;
          }
          pff = pff->next;
          k++;
        }
      }
    ff = ff->next;
    gauge.GaugeStep( l++, 0 );
  }
  gauge.CloseGauge();
}
//---------------------------------------------------------------------------
CovMatrix * CreateMatrix( TBoolFnct * fn )
{
  CovMatrix * cm;
  int i, j;
  TPfnct * ff;
  TCovs * c;
  CovLink * cl;

  cm = new CovMatrix;
  cm->matrix = new char * [fn->terms];
  for ( i = 0; i < fn->terms; i++ ) cm->matrix[i] = new char[fn->parent->maxonterms * fn->outvars];
  cm->covs = 0;
  cm->cols = fn->parent->maxonterms * fn->outvars;
  cm->rows = fn->terms;
  cm->origrows = cm->rows;
  cm->solcost = 0;
  cm->col = 0;
  cm->row = 0;
  cm->sol = 0;
  cm->cleancheck = 0;
  cm->fn = fn;
  cm->sparsing = 0;

  for ( i = 0; i < fn->outvars; i++ )
  {
    ff = fn->parent->onset[i];
    for ( j = 0; j < fn->parent->onterms[i]; j++ )
    {
      cl = new CovLink;
      if ( cm->col != 0 ) cm->col->prev = cl;
      cl->index = fn->parent->maxonterms * i + j;
      cl->term = ff;
      cl->next = cm->col;
      cl->prev = 0;
      cm->col = cl;
      ff = ff->next;
    }
    cm->covs += j;
  }
  ff = fn->fnct;
  i = 0;
  while ( ff != 0 )
  {
    cl = new CovLink;
    if ( cm->row != 0 ) cm->row->prev = cl;
    cl->index = i;
    cl->term = ff;
    cl->next = cm->row;
    cl->prev = 0;
    cm->row = cl;
    for ( j = 0; j < fn->parent->maxonterms * fn->outvars; j++ ) cm->matrix[i] [j] = '0';
    c = ff->cov;
    while ( c != 0 )
    {
      cm->matrix[i] [c->var * fn->parent->maxonterms + c->index] = '1';
      c = c->next;
    }
    i++;
    ff = ff->next;
  }
  return cm;
}
//---------------------------------------------------------------------------
CovMatrix * CreateBlockMatrix( TBoolFnct * fn, int block )
{
  CovMatrix * cm;
  int i, j;
  TPfnct * ff;
  TCovs * c;
  CovLink * cl;

  cm = new CovMatrix;
  cm->matrix = new char * [fn->terms];
  for ( i = 0; i < fn->terms; i++ ) cm->matrix[i] = new char[fn->parent->maxonterms * fn->outvars];
  cm->covs = 0;
  cm->cols = fn->parent->maxonterms * fn->outvars;
  cm->rows = fn->terms;
  cm->origrows = cm->rows;
  cm->solcost = 0;
  cm->col = 0;
  cm->row = 0;
  cm->sol = 0;
  cm->cleancheck = 0;
  cm->fn = fn;
  cm->block = block;
  cm->sparsing = 0;

  for ( i = 0; i < fn->outvars; i++ )
    if ( fn->decomp->outs[i] == block )
    {
      ff = fn->parent->onset[i];
      for ( j = 0; j < fn->parent->onterms[i]; j++ )
      {
        cl = new CovLink;
        if ( cm->col != 0 ) cm->col->prev = cl;
        cl->index = fn->parent->maxonterms * i + j;
        cl->term = ff;
        cl->next = cm->col;
        cl->prev = 0;
        cm->col = cl;
        ff = ff->next;
      }
      cm->covs += j;
    }
  ff = fn->fnct;
  i = 0;
  while ( ff != 0 )
  {
    cl = new CovLink;
    if ( cm->row != 0 ) cm->row->prev = cl;
    cl->index = i;
    cl->term = ff;
    cl->next = cm->row;
    cl->prev = 0;
    cm->row = cl;
    for ( j = 0; j < fn->parent->maxonterms * fn->outvars; j++ ) cm->matrix[i] [j] = '0';
    c = ff->cov;
    while ( c != 0 )
    {
      if ( fn->decomp->outs[c->var] == block )
        cm->matrix[i] [c->var * fn->parent->maxonterms + c->index] = '1';
      c = c->next;
    }
    i++;
    ff = ff->next;
  }
  return cm;
}
//---------------------------------------------------------------------------
CovMatrix * CreateSOMatrix( TBoolFnct * fn, int ovr )
{
  CovMatrix * cm;
  int i, j;
  TPfnct * ff;
  TCovs * c;
  CovLink * cl;

  cm = new CovMatrix;
  cm->matrix = new char * [fn->terms];
  for ( i = 0; i < fn->terms; i++ ) cm->matrix[i] = new char[fn->parent->onterms[ovr]];
  cm->cols = fn->parent->onterms[ovr];
  cm->rows = fn->terms;
  cm->origrows = cm->rows;
  cm->solcost = 0;
  cm->col = 0;
  cm->row = 0;
  cm->sol = 0;
  cm->cleancheck = 0;
  cm->fn = fn;
  cm->sparsing = 0;

  ff = fn->parent->onset[ovr];
  for ( j = 0; j < fn->parent->onterms[ovr]; j++ )
  {
    cl = new CovLink;
    if ( cm->col != 0 ) cm->col->prev = cl;
    cl->index = j;
    cl->term = ff;
    cl->next = cm->col;
    cl->prev = 0;
    cm->col = cl;
    ff = ff->next;
  }
  cm->covs = j;

  ff = fn->fnct;
  i = 0;
  while ( ff != 0 )
  {
    cl = new CovLink;
    if ( cm->row != 0 ) cm->row->prev = cl;
    cl->index = i;
    cl->term = ff;
    cl->next = cm->row;
    cl->prev = 0;
    cm->row = cl;
    for ( j = 0; j < fn->parent->onterms[ovr]; j++ ) cm->matrix[i] [j] = '0';
    c = ff->cov;
    while ( c != 0 )
    {
      if ( c->var == ovr && ff->ov[ovr] == '1' ) cm->matrix[i] [c->index] = '1';
      c = c->next;
    }
    i++;
    ff = ff->next;
  }
  return cm;
}
//---------------------------------------------------------------------------
void DelLinks( CovLink * cl )
{
  CovLink * cl1;

  while ( cl != 0 )
  {
    cl1 = cl->next;
    delete cl;
    cl = cl1;
  }
}
//---------------------------------------------------------------------------
void DeleteMatrix( CovMatrix * cm )
{
  int i;

  DelLinks( cm->row );
  DelLinks( cm->col );
  for ( i = 0; i < cm->origrows; i++ ) delete[] cm->matrix[i];
  delete[] cm->matrix;
  delete cm;
}
//---------------------------------------------------------------------------
CovMatrix * ForkMatrix( CovMatrix * cm )
{
  CovMatrix * cm1;
  CovLink * cl, * cl1;

  cm1 = new( CovMatrix );
  cm1->matrix = cm->matrix;
  cm1->col = 0;
  cm1->row = 0;
  cm1->sol = 0;
  cm1->cleancheck = 0;

  cl = cm->col;
  while ( cl != 0 )
  {
    cl1 = new( CovLink );
    if ( cm1->col != 0 ) cm1->col->prev = cl1;
    cl1->index = cl->index;
    cl1->next = cm1->col;
    cl1->prev = 0;
    cm1->col = cl1;
    cl = cl->next;
  }
  cl = cm->row;
  while ( cl != 0 )
  {
    cl1 = new( CovLink );
    if ( cm1->row != 0 ) cm1->row->prev = cl1;
    cl1->index = cl->index;
    cl1->next = cm1->row;
    cl1->prev = 0;
    cm1->row = cl1;
    cl = cl->next;
  }
  cl = cm->sol;
  while ( cl != 0 )
  {
    cl1 = new( CovLink );
    if ( cm1->sol != 0 ) cm1->sol->prev = cl1;
    cl1->index = cl->index;
    cl1->next = cm1->sol;
    cl1->prev = 0;
    cm1->sol = cl1;
    cl = cl->next;
  }
  cm1->covs = cm->covs;
  cm1->cols = cm->cols;
  cm1->rows = cm->rows;
  cm1->solcost = cm->solcost;
  return cm1;
}
//---------------------------------------------------------------------------
void SelRow( CovMatrix * cm, CovLink * r )
{
  CovLink * cl, * cl1, * cl2;

  if ( r->prev == 0 )
  {
    cm->row = r->next;
    if ( r->next != 0 ) r->next->prev = 0;
  }
  else
  {
    r->prev->next = r->next;
    if ( r->next != 0 ) r->next->prev = r->prev;
  }
  r->prev = 0;
  if ( cm->sol != 0 ) cm->sol->prev = r;
  r->next = cm->sol;
  cm->sol = r;
  cm->solcost++;

  cl = cm->col;
  while ( cl != 0 )
  {
    if ( cm->matrix[r->index] [cl->index] == '1' )
    {
      cm->covs--;
      cl1 = cl->next;
      if ( cl->prev == 0 )
      {
        cm->col = cl->next;
        if ( cl->next != 0 ) cl->next->prev = 0;
      }
      else
      {
        cl->prev->next = cl->next;
        if ( cl->next != 0 ) cl->next->prev = cl->prev;
      }
      if ( cm->cleancheck != 0 )
      {
        cl2 = cm->row;
        while ( cl2 != 0 )
        {
          if ( cm->matrix[cl2->index] [cl->index] == '1' )
          {
            if ( cm->cleancheck == 2 ) cl2->contrib -= 1 / cl->contrib;
            //cl2->contrib-=1/((cl->contrib!=1)?cl->contrib-1:0.00001);
            //cl2->contrib-=(cm->rows-cl->contrib+0.00001)/((cl->contrib!=1)?cl->contrib-1:0.00001);
            cl2->count--;
          }
          cl2 = cl2->next;
        }
      }
      cm->cols--;
      delete cl;
      cl = cl1;
    }
    else
      cl = cl->next;
  }
  cm->rows--; // to je tady novy!
}
//---------------------------------------------------------------------------
void CutRow( CovMatrix * cm, CovMatrix * cm1, CovLink * r )
{
  CovLink * cl, * cl1;

  if ( r->prev == 0 )
  {
    cm->row = r->next;
    if ( r->next != 0 ) r->next->prev = 0;
  }
  else
  {
    r->prev->next = r->next;
    if ( r->next != 0 ) r->next->prev = r->prev;
  };
  r->prev = 0;
  if ( cm1->row != 0 ) cm1->row->prev = r;
  r->next = cm1->row;
  cm1->row = r;

  cl = cm->col;
  while ( cl != 0 )
  {
    if ( cm->matrix[r->index] [cl->index] == '1' )
    {
      cm->covs--;
      cm1->covs++;
      cl1 = cl->next;
      if ( cl->prev == 0 )
      {
        cm->col = cl->next;
        if ( cl->next != 0 ) cl->next->prev = 0;
      }
      else
      {
        cl->prev->next = cl->next;
        if ( cl->next != 0 ) cl->next->prev = cl->prev;
      };
      cl->prev = 0;
      if ( cm1->col != 0 ) cm1->col->prev = cl;
      cl->next = cm1->col;
      cm1->col = cl;
      cl = cl1;
    }
    else
      cl = cl->next;
  }
}
//---------------------------------------------------------------------------
void DelRow( CovMatrix * cm, CovLink * r )
{
  CovLink * cl, * cl1;

  if ( r->prev == 0 )
  {
    cm->row = r->next;
    if ( r->next != 0 ) r->next->prev = 0;
  }
  else
  {
    r->prev->next = r->next;
    if ( r->next != 0 ) r->next->prev = r->prev;
  };
  if ( cm->cleancheck != 0 )
  {
    cl = cm->col;
    while ( cl != 0 )
    {
      if ( cm->matrix[r->index] [cl->index] == '1' )
      {
        if ( cm->cleancheck == 2 )
        {
          cl1 = cm->row;
          while ( cl1 != 0 )
          {
            if ( cm->matrix[cl1->index] [cl->index] == '1' )
            {
              cl1->contrib -= 1 / cl->contrib;
              cl1->contrib += 1 / ( cl->contrib - 1 );
            }
            cl1 = cl1->next;
          }
        }
        cl->contrib--;
        cl->count--;
      }
      cl = cl->next;
    }
  }
  delete r;
}
//---------------------------------------------------------------------------
CovMatrix * SplitMatrix( CovMatrix * cm, char * cycl )
{
  CovLink * cl;
  CovMatrix * cm1;

  cm1 = new( CovMatrix );
  cm1->matrix = cm->matrix;
  cm1->col = 0;
  cm1->row = 0;
  cm1->sol = 0;
  cm1->cleancheck = 0;

  cm1->cols = cm->cols;
  cm1->rows = cm->rows;
  cm1->solcost = 0;
  cm1->covs = 0;

  cl = cm->row;
  while ( cl != 0 )
  {
    if ( cycl[cl->index] == '1' )
    {
      CutRow( cm, cm1, cl );
      cl = cm->row;
    }
    else
      cl = cl->next;
  }
  return cm1;
}

void SaveMatrix( CovMatrix * cm )
{
  CovLink * cl, * cl1;
  FILE * f;
  int i, j;

  f = fopen( "matrix.dat", "w" );
  cl = cm->row;
  while ( cl != 0 )
  {
    fprintf( f, "%i:", cl->index );
    cl1 = cm->col;
    while ( cl1 != 0 )
    {
      fprintf( f, "%i ", cl1->index );
      cl1 = cl1->next;
    }
    fprintf( f, "\n" );
    cl = cl->next;
  }
  fprintf( f, "\n" );

  for ( i = 0; i < cm->rows; i++ )
  {
    for ( j = 0; j < cm->cols; j++ ) fprintf( f, "%c", cm->matrix[i] [j] );
    fprintf( f, "\n" );
  }
  fprintf( f, "\n" );

  cl = cm->row;
  while ( cl != 0 )
  {
    cl1 = cm->col;
    while ( cl1 != 0 )
    {
      fprintf( f, "%c", cm->matrix[cl->index] [cl1->index] );
      cl1 = cl1->next;
    }
    fprintf( f, "\n" );
    cl = cl->next;
  }
  fprintf( f, "\n" );

  fclose( f );
}

void SaveMatrix2(CovMatrix *cm)
{
    CovLink *cl, *cl1;
    FILE *f;
    int i,j;

    f=fopen("matrix.dat","w");
/*
    cl=cm->row;
    while (cl!=NULL) {
        fprintf(f,"%i:",cl->index);
        cl1=cm->col;
        while (cl1!=NULL) {
            fprintf(f,"%i ",cl1->index);
            cl1=cl1->next;
        }
        fprintf(f,"\n");
        cl=cl->next;
    }
    fprintf(f,"\n");

    for (i=0; i<cm->rows; i++) {
        for (j=0; j<cm->cols; j++) fprintf(f,"%c",cm->matrix[i][j]);
        fprintf(f,"\n");
    }
    fprintf(f,"\n");
*/
    fprintf(f, "%i %i\n", cm->rows, cm->cols);
    cl1 = cm->row;
    while (cl1 !=NULL ) {
        fprintf(f, "%i ", cm->fn->inpvars - cl1->term->dim);
        cl1 = cl1->next;
    }
    fprintf(f, "\n");

    cl=cm->col;
    while (cl!=NULL) {
        cl1=cm->row;
        while (cl1!=NULL) {
            fprintf(f,"%c ",cm->matrix[cl1->index][cl->index]);
            cl1=cl1->next;
        }
        fprintf(f,"\n");
        cl=cl->next;
    }
    fprintf(f,"\n");

    fclose(f);
}

//---------------------------------------------------------------------------
void ColDomChck( CovMatrix * cm )
{
  int flg;
  CovLink * cl, * cl1, * r;

  cl = cm->col;
  while ( cl != 0 )
  {
    flg = 0;
    cl1 = cm->row;
    while ( cl1 != 0 )
    {
      if ( cm->matrix[cl1->index] [cl->index] == '1' )
      {
        flg++;
        r = cl1;
      }
      if ( flg > 1 ) break;
      cl1 = cl1->next;
    }
    if ( flg == 1 )
    {
      SelRow( cm, r );
      cl = cm->col;
    }
    else
      cl = cl->next;
  }
}
//---------------------------------------------------------------------------
void RowDomChck( CovMatrix * cm )
{
  int flg, dflg;
  CovLink * cl, * cl1, * cl2, * tcl;

  cl = cm->row;
  while ( cl != 0 )
  {
    cl1 = cm->row;
    dflg = 0;
    while ( cl1 != 0 )
    {
      if ( cl1 != cl )
      {
        flg = 0;
        cl2 = cm->col;
        while ( cl2 != 0 )
        {
          if ( cm->matrix[cl->index] [cl2->index] == '1' && cm->matrix[cl1->index] [cl2->index] == '0' ) flg = 1;
          if ( flg == 1 ) break;
          cl2 = cl2->next;
        };
        if ( flg == 0 )
        {
          tcl = cl->next;
          DelRow( cm, cl );
          cl = tcl;
          dflg = 1;
          cm->dact = 1;
          break;
        }
      }
      cl1 = cl1->next;
    }
    if ( dflg == 0 ) cl = cl->next;
  }
}
//---------------------------------------------------------------------------
void CreateVct( CovMatrix * cm )
{
  CovLink * cl, * cl1;
  int p;
  unsigned long b;

  cl = cm->row;
  while ( cl != 0 )
  {
    memset( cl->vct, 0, cm->rows / sizeof( unsigned long ) + 1 );
    cl1 = cm->col;
    p = 0;
    b = 1;
    while ( cl1 != 0 )
    {
      if ( cm->matrix[cl->index] [cl1->index] == '1' ) cl->vct[p] |= b;
      b <<= 1;
      if ( b == 0 )
      {
        b = 1;
        p++;
      }
      cl1 = cl1->next;
    }
    cl = cl->next;
  }
}
//---------------------------------------------------------------------------
void RowDomChck2( CovMatrix * cm )
{
  int dflg;
  CovLink * cl, * cl1, * tcl;
  unsigned int p;
  unsigned long r1, r2, tr;

  CreateVct( cm );
  cl = cm->row;
  while ( cl != 0 )
  {
    cl1 = cm->row;
    dflg = 0;
    while ( cl1 != 0 )
    {
      if ( cl1 != cl )
      {
        r1 = 0;
        r2 = 0;
        for ( p = 0; p < cm->rows / sizeof( unsigned long ) + 1; p++ )
        {
          tr = cl->vct[p] ^ cl1->vct[p];
          r1 |= tr & cl->vct[p];
          r2 |= tr & cl1->vct[p];
        };
        if ( r1 == 0 && r2 != 0 )
        {
          SaveMatrix( cm );
          tcl = cl->next;
          DelRow( cm, cl );
          cl = tcl;
          dflg = 1;
          cm->dact = 1;
          CreateVct( cm );
          break;
        }
      }
      cl1 = cl1->next;
    }
    if ( dflg == 0 ) cl = cl->next;
  }
}
//---------------------------------------------------------------------------
int MSIC( CovMatrix * cm, char * crows, CovLink * pos, int m )
{
  int j;
  char * crows1;
  int m1, m2;
  CovLink * cl;

  while ( pos != 0 )
  {
    cl = cm->row;
    while ( cl != 0 )
    {
      if ( cm->matrix[cl->index] [pos->index] == '1' && crows[cl->index] == '1' ) break;
      cl = cl->next;
    }
    if ( cl == 0 ) break;
    pos = pos->next;
  }
  if ( pos == 0 ) return m;
  m1 = max( m, MSIC( cm, crows, pos->next, m ) );
  crows1 = new( char[cm->rows] );
  for ( j = 0; j < cm->rows; j++ ) crows1[j] = crows[j];
  for ( j = 0; j < cm->rows; j++ )
    if ( cm->matrix[j] [pos->index] == '1' ) crows1[j] = '1';
  m2 = max( m1, MSIC( cm, crows1, pos->next, m + 1 ) );
  delete[] crows1;
  return m2;
}
//---------------------------------------------------------------------------
int FindMSIC( CovMatrix * cm )
{
  char * crows;
  int i;
  int m;

  crows = new( char[cm->rows] );
  for ( i = 0; i < cm->rows; i++ ) crows[i] = '0';

  m = MSIC( cm, crows, cm->col, 0 );
  delete[] crows;
  return m;
}
//---------------------------------------------------------------------------
int HeurFindMSIC2( CovMatrix * cm )
{
  char * crows;
  int i, j;
  int m;
  CovLink * cl;
  CovLink * pos;

  crows = new( char[cm->rows] );
  for ( i = 0; i < cm->rows; i++ ) crows[i] = '0';
  m = 0;
  pos = cm->col;

  do
  {
    while ( pos != 0 )
    {
      cl = cm->row;
      while ( cl != 0 )
      {
        if ( cm->matrix[cl->index] [pos->index] == '1' && crows[cl->index] == '1' ) break;
        cl = cl->next;
      }
      if ( cl == 0 ) break;
      pos = pos->next;
    }
    if ( pos == 0 )
    {
      delete[] crows;
      return m;
    }
    m++;
    for ( j = 0; j < cm->rows; j++ )
      if ( cm->matrix[j] [pos->index] == '1' ) crows[j] = '1';
    pos = pos->next;
  }
  while ( 1 );
}
//---------------------------------------------------------------------------
int HeurFindMSIC( CovMatrix * cm )
{
  char * crows;
  int i, j;
  int m;
  CovLink * cl, * cl1;
  CovLink * pos;
  int min;

  crows = new( char[cm->rows] );
  for ( i = 0; i < cm->rows; i++ ) crows[i] = '0';
  m = 0;
  cl = cm->col;
  while ( cl != 0 )
  {
    cl->count = 0;
    cl1 = cm->row;
    while ( cl1 != 0 )
    {
      if ( cm->matrix[cl1->index] [cl->index] == '1' ) cl->count++;
      cl1 = cl1->next;
    };
    cl = cl->next;
  }
  do
  {
    pos = cm->col;
    min = INT_MAX;
    cl1 = 0;
    while ( pos != 0 )
    {
      cl = cm->row;
      while ( cl != 0 )
      {
        if ( cm->matrix[cl->index] [pos->index] == '1' && crows[cl->index] == '1' ) break;
        cl = cl->next;
      }
      if ( cl == 0 && pos->count < min )
      {
        min = pos->count;
        cl1 = pos;
      }
      pos = pos->next;
    }
    if ( cl1 == 0 )
    {
      delete[] crows;
      return m;
    }
    pos = cl1;
    m++;
    for ( j = 0; j < cm->rows; j++ )
      if ( cm->matrix[j] [pos->index] == '1' ) crows[j] = '1';
    pos = pos->next;
  }
  while ( 1 );
}
//---------------------------------------------------------------------------
int F_Cycl( CovMatrix * cm, char * cycl, CovLink * pos, int num )
{
  CovLink * cl, * cl1;

  cycl[pos->index] = '1';
  num++;
  cl = cm->col;
  while ( cl != 0 )
  {
    if ( cm->matrix[pos->index] [cl->index] == '1' )
    {
      cl1 = cm->row;
      while ( cl1 != 0 )
      {
        if ( cm->matrix[cl1->index] [cl->index] == '1' && cycl[cl1->index] == '0' )
          num = F_Cycl( cm, cycl, cl1, num );
        cl1 = cl1->next;
      }
    }
    cl = cl->next;
  }
  return num;
}
//---------------------------------------------------------------------------
char * FindCycle( CovMatrix * cm )
{
  char * cycl;
  int i, j;
  int c;
  CovLink * cl;

  cycl = new( char[cm->rows] );
  for ( i = 0; i < cm->rows; i++ ) cycl[i] = '0';

  c = F_Cycl( cm, cycl, cm->row, 0 );

  j = 0;
  cl = cm->row;
  while ( cl != 0 )
  {
    j++;
    cl = cl->next;
  }
  if ( c == j ) return 0;
  return cycl;
}
//---------------------------------------------------------------------------
void CovBranch( CovMatrix * cm, CovSol * best )
{
  CovMatrix * cm1;
  int dcovs;
  int i;
  CovLink * cl;
  char * cycl;
  CovSol * best1, * best2;

  cm->cleancheck = 0;

  if ( cm->solcost >= best->solcost )
  {
    DelLinks( cm->sol );
    return;
  }
  do
  {
    dcovs = cm->covs;
    ColDomChck( cm );
    dcovs -= cm->covs;
    if ( cm->solcost >= best->solcost )
    {
      DelLinks( cm->sol );
      return;
    }
    cm->dact = 0;
    if ( cm->covs > 0 ) RowDomChck( cm );
  }
  while ( cm->covs > 0 && ( dcovs > 0 || cm->dact == 1 ) );
  if ( cm->covs == 0 )
  {
    if ( cm->solcost < best->solcost )
    {
      DelLinks( best->sol );
      best->sol = cm->sol;
      best->solcost = cm->solcost;
    }
    else
      DelLinks( cm->sol );
    return;
  }
  if ( HeurFindMSIC( cm ) + cm->solcost >= best->solcost ) return;
  cycl = FindCycle( cm );
  if ( cycl != 0 )
  {
    cm1 = SplitMatrix( cm, cycl );
    best1 = new( CovSol );
    best1->solcost = INT_MAX;
    best1->sol = 0;
    CovBranch( cm, best1 );
    best2 = new( CovSol );
    best2->solcost = INT_MAX;
    best2->sol = 0;
    CovBranch( cm1, best2 );
    best1->solcost += best2->solcost;
    if ( best1->solcost < best->solcost )
    {
      DelLinks( best->sol );
      cl = best2->sol;
      while ( cl->next != 0 ) cl = cl->next;
      best->sol = best2->sol;
      cl->next = best1->sol;
      best1->sol->prev = cl;
    }
    else
    {
      DelLinks( best1->sol );
      DelLinks( best2->sol );
    }
    delete best1;
    delete best2;
  }
  else
  {
    cm1 = ForkMatrix( cm );
    i = cm1->row->index;
    SelRow( cm1, cm1->row );
    CovBranch( cm1, best );
    DelLinks( cm1->col );
    DelLinks( cm1->row );
    delete cm1;
    cl = cm->row;
    while ( cl->index != i ) cl = cl->next;
    DelRow( cm, cl );
    CovBranch( cm, best );
  }
  return;
}
//---------------------------------------------------------------------------
void PriorCover( CovMatrix * cm, CovSol * best, TGauge gauge )
{
  int j, prg;
  CovLink * cl, * cl1, * scl;
  double min, max, mdim;
  int minc, maxc, mdimc;

  gauge.NewGauge( "Solving CP - Priority cover", cm->covs, 1 );
  prg = cm->covs;
  cm->cleancheck = 1;
  cl = cm->row;
  while ( cl != 0 )
  {
    cl->count = 0;
    cl = cl->next;
  }
  cl = cm->col;
  while ( cl != 0 )
  {
    cl->count = 0;
    cl = cl->next;
  }
  cl = cm->col; // pocitani cetnosti
  while ( cl != 0 )
  {
    cl1 = cm->row;
    while ( cl1 != 0 )
    {
      if ( cm->matrix[cl1->index] [cl->index] == '1' )
      {
        cl1->count++;
        cl->count++;
      }
      cl1 = cl1->next;
    }
    cl = cl->next;
  }
  do
  {
    //RowDomChck2(cm);

    min = INT_MAX; // minimalni pocet radku
    minc = 0;
    cl1 = cm->col;
    while ( cl1 != 0 )
    {
      if ( cl1->count == min ) minc++;
      if ( cl1->count < min )
      {
        min = cl1->count;
        minc = 1;
      }
      cl1 = cl1->next;
    }
    max = 0; // maximalni pocet sloupcu
    maxc = 0;
    cl = cm->col;
    while ( cl != 0 )
    {
      if ( cl->count == min )
      {
        cl1 = cm->row;
        while ( cl1 != 0 )
        {
          if ( cm->matrix[cl1->index] [cl->index] == '1' )
          {
            if ( cl1->count == max ) maxc++;
            if ( cl1->count > max )
            {
              max = cl1->count;
              maxc = 1;
              scl = cl1;
            }
          }
          cl1 = cl1->next;
        }
      }
      cl = cl->next;
    }
    if ( maxc != 1 )
    {
      cl = cm->col; // nejvetsi dimenze
      mdim = INT_MAX;
      mdimc = 0;
      while ( cl != 0 )
      {
        if ( cl->count == min )
        {
          cl1 = cm->row;
          while ( cl1 != 0 )
          {
            if ( cm->matrix[cl1->index] [cl->index] == '1' && cl1->count == max )
            {
              if ( cl1->term->dim == mdim ) mdimc++;
              if ( cl1->term->dim > mdim )
              {
                mdim = cl1->term->dim;
                mdimc = 1;
                scl = cl1;
              }
            }
            cl1 = cl1->next;
          }
        }
        cl = cl->next;
      }
    }
    if ( mdimc != 1 )
    {
      cl = cm->col; // nahodny vyber
      j = random( mdimc ) + 1;
      while ( cl != 0 )
      {
        if ( cl->count == min )
        {
          cl1 = cm->row;
          while ( cl1 != 0 )
          {
            if ( cm->matrix[cl1->index] [cl->index] == '1' && cl1->count == max && cl1->term->dim == mdim )
            {
              j--;
              scl = cl1;
              if ( j == 0 ) break;
            }
            cl1 = cl1->next;
          }
        }
        if ( j == 0 ) break;
        cl = cl->next;
      }
    }
    SelRow( cm, scl );
    gauge.GaugeStep( prg - cm->covs, 0 );
  }
  while ( cm->covs != 0 );
  best->sol = cm->sol;
  best->solcost = cm->solcost;
  gauge.CloseGauge();
}
//---------------------------------------------------------------------------
void CountContrib( CovMatrix * cm )
{
  CovLink * cl, * cl1;

  cl1 = cm->row;
  while ( cl1 != 0 )
  {
    cl1->contrib = 0;
    cl1 = cl1->next;
  }
  cl = cm->col;
  while ( cl != 0 )
  {
    cl1 = cm->row;
    cl->contrib = 0;
    while ( cl1 != 0 )
    {
      if ( cm->matrix[cl1->index] [cl->index] == '1' ) cl->contrib++;
      cl1 = cl1->next;
    }
    cl1 = cm->row;
    while ( cl1 != 0 )
    {
      if ( cm->matrix[cl1->index] [cl->index] == '1' ) cl1->contrib += 1 / cl->contrib;
      //cl1->contrib+=1/((cl->contrib!=1)?cl->contrib-1:0.00001);
      //cl1->contrib+=(cm->rows-cl->contrib+0.00001)/((cl->contrib!=1)?cl->contrib-1:0.00001);
      cl1 = cl1->next;
    }
    cl = cl->next;
  }
}
//---------------------------------------------------------------------------
void ContribCover( CovMatrix * cm, CovSol * best, TGauge gauge )
{
  CovLink * cl1, * scl;
  double min;
  int prg;

  gauge.NewGauge( "Solving CP - COntributive removal", cm->covs, 1 );
  cm->cleancheck = 2;
  CountContrib( cm );
  prg = cm->covs;
  do
  {
    //RowDomChck(cm);
    ColDomChck( cm );
    if ( cm->covs != 0 )
    {
      cl1 = cm->row;
      min = LONG_MAX;
      while ( cl1 != 0 )
      {
        if ( cl1->contrib < min )
        {
          min = cl1->contrib;
          scl = cl1;
        }
        cl1 = cl1->next;
      }
      DelRow( cm, scl );
    }
    gauge.GaugeStep( prg - cm->covs, 0 );
  }
  while ( cm->covs != 0 );
  best->sol = cm->sol;
  best->solcost = cm->solcost;
  gauge.CloseGauge();
}
//---------------------------------------------------------------------------
void ContribAddCover( CovMatrix * cm, CovSol * best, TGauge gauge, Cov_Opts * covopts )
{
  int i, j, cnt;
  CovLink * cl1, * scl;
  double wcontrib, max;
  char * * iusg;
  int prg;

  gauge.NewGauge( "Solving CP - Contributive addition", cm->covs, 1 );
  prg = cm->covs;

  if ( covopts->DFT_COVenable == 1 && cm->sparsing == 0 )
  {
    iusg = new char * [cm->fn->outvars];
    for ( i = 0; i < cm->fn->outvars; i++ )
    {
      iusg[i] = new char[cm->fn->inpvars];
      setmem( iusg[i], cm->fn->inpvars, 0 );
    }
  }

  cm->cleancheck = 2;
  CountContrib( cm );
  while ( cm->covs != 0 )
  {
    cl1 = cm->row;
    max = 0;
    while ( cl1 != 0 )
    {
      if ( cl1->contrib > 0 )
      {
        cnt = 0;
        if ( covopts->DECOMP_COVenable == 1 && cm->sparsing == 0 )
        {
          for ( i = 0; i < cm->fn->inpvars; i++ )
            if ( cl1->term->iv[i] != '-' && cm->fn->decomp->inps[cm->block] [i] == 'N' ) cnt++;
        }
        if ( covopts->DECOMP_COVenable == 1 ) wcontrib = cl1->contrib / ( ( ( double )cnt * covopts->DECOMP_COV_force ) + 1 );
        else
          wcontrib = cl1->contrib;
        if ( covopts->DECOMP_COVenable == 1 && covopts->DECOMP_COV_loadforce > 0 )
        {
          cnt = 0;
          for ( j = 0; j < cm->fn->inpvars; j++ )
            if ( cl1->term->iv[j] != '-' )
              for ( i = 0; i < cm->fn->decomp->blocks; i++ )
                if ( i != cm->block && cm->fn->decomp->inps[i] [j] != 'N' )
                {
                  cnt++;
                  break;
                }
          wcontrib /= ( ( ( double )cnt * covopts->DECOMP_COV_loadforce ) + 1 );
        }
        if ( covopts->DFT_COVenable == 1 && cm->sparsing == 0 )
        {
          cnt = 0;
          for ( j = 0; j < cm->fn->outvars; j++ )
            for ( i = 0; i < cm->fn->inpvars; i++ )
              if ( cl1->term->ov[j] == '1' && cl1->term->iv[i] != '-' && iusg[j] [i] == 0 ) cnt++;
          wcontrib = wcontrib / ( ( ( double )cnt * covopts->DFT_COV_force ) + 1 );
        }
        if ( wcontrib > max )
        {
          max = wcontrib;
          scl = cl1;
        }
      }
      cl1 = cl1->next;
    }
    SelRow( cm, scl );
    if ( covopts->DECOMP_COVenable == 1 && cm->sparsing == 0 )
      for ( i = 0; i < cm->fn->inpvars; i++ )
        if ( scl->term->iv[i] != '-' )
        {
          if ( cm->fn->decomp->inps[cm->block] [i] == 'N' ) cm->fn->decomp->inps[cm->block] [i] = '1';
          else if ( cm->fn->decomp->inps[cm->block] [i] == 'R' ) cm->fn->decomp->inps[cm->block] [i] = 'F';
        }
    if ( covopts->DFT_COVenable == 1 && cm->sparsing == 0 )
    {
      for ( j = 0; j < cm->fn->outvars; j++ )
        for ( i = 0; i < cm->fn->inpvars; i++ )
          if ( scl->term->ov[j] == '1' && scl->term->iv[i] != '-' ) iusg[j] [i] = 1;
    }
    gauge.GaugeStep( prg - cm->covs, 0 );
  }
  best->sol = cm->sol;
  best->solcost = cm->solcost;

  if ( covopts->DFT_COVenable == 1 && cm->sparsing == 0 )
  {
    for ( i = 0; i < cm->fn->outvars; i++ )
      delete[] iusg[i];
    delete[] iusg;
  }
  gauge.CloseGauge();
}
//---------------------------------------------------------------------------
CovSol * ProcessCover( CovMatrix * cm, Cov_Opts * covopts, TGauge gauge )
{
  CovSol * best;

  best = new CovSol;
  best->solcost = INT_MAX;
  best->sol = 0;
  if ( cm->cols == 0 ) return best;
  switch ( covopts->COV_type )
  {
    case 0:
      PriorCover( cm, best, gauge);
    break;
    case 1:
      ContribAddCover( cm, best, gauge, covopts );
    break;
    case 2:
      ContribCover( cm, best, gauge );
    break;
    case 3:
      gauge.NewGauge( "Solving CP - exact", 0, 1 );
      CovBranch( cm, best );
      gauge.CloseGauge();
    break;
  }
  return best;
}
//---------------------------------------------------------------------------
TBoolFnct * Cover( TBoolFnct * fn, TGauge gauge, Cov_Opts * covopts )
{
  CovMatrix * cm;
  TBoolFnct * outfn;
  TPfnct * ff, * ff1;
  CovLink * cl;
  CovSol * best;

  cm = CreateMatrix( fn );
//  SaveMatrix2(cm);
  best = ProcessCover( cm, covopts, gauge );
  outfn = fn->parent->CreateChildFnct();
  outfn->type = 0;

  cl = best->sol;
  while ( cl != 0 )
  {
    ff = cl->term;
    ff1 = outfn->AddTermToEnd( ff->iv, ff->ov, ff->dim );
    ff1->cov = DuplCov( ff->cov );
    ff1->covlen = ff->covlen;
    if ( ff->covfield != 0 )
    {
      ff1->covfield = new int[ff->covlen];
      memcpy( ff1->covfield, ff->covfield, ff->covlen * sizeof( int ) );
    }
    ff1->prime = ff->prime;
    ff1->origin = ff->origin;
    cl = cl->next;
  }
  DelLinks( best->sol );
  delete best;
  DeleteMatrix( cm );
  return outfn;
}
//---------------------------------------------------------------------------
TBoolFnct * SOCover( TBoolFnct * fn, Cov_Opts * covopts, TGauge gauge )
{
  int i;
  CovMatrix * cm;
  TBoolFnct * outfn;
  TPfnct * ff1;
  CovLink * cl;
  CovSol * best;
  char * ov;

  outfn = fn->parent->CreateChildFnct();
  outfn->type = 0;
  ov = new char[fn->outvars + 1];
  for ( i = 0; i < fn->outvars; i++ ) ov[i] = '0';
  ov[i] = 0;

  for ( i = 0; i < fn->outvars; i++ )
  {
    cm = CreateSOMatrix( fn, i );
    best = ProcessCover( cm, covopts, gauge );
    cl = best->sol;
    while ( cl != 0 )
    {
      ff1 = outfn->AddTermToEnd( cl->term->iv, ov, cl->term->dim );
      ff1->ov[i] = '1';
      cl = cl->next;
    }
    DelLinks( best->sol );
    delete best;
    DeleteMatrix( cm );
  }
  return outfn;
}
//---------------------------------------------------------------------------
void SparseOutput( TBoolFnct * fn, Cov_Opts * covopts, TGauge gauge )
{
  int i;
  CovMatrix * cm;
  TPfnct * ff;
  CovLink * cl;
  CovSol * best;

  gauge.NewGauge( "Sparsing output matrix", fn->outvars, 1 );
  for ( i = 0; i < fn->outvars; i++ )
  {
    cm = CreateSOMatrix( fn, i );
    cm->sparsing = 1;
    best = ProcessCover( cm, covopts, gauge );
    ff = fn->fnct;
    while ( ff != 0 )
    {
      if ( ff->ov[i] == '1' )
      {
        cl = best->sol;
        while ( cl != 0 )
        {
          if ( cl->term == ff ) break;
          cl = cl->next;
        }
        if ( cl == 0 )
        {
          if ( covopts->COV_substDC == 1 ) ff->ov[i] = '0';
          else
            ff->ov[i] = '~';
        }
      }
      ff = ff->next;
    }
    DelLinks( best->sol );
    delete best;
    DeleteMatrix( cm );
    gauge.GaugeStep( i, 0 );
  }
  // delete useless terms
  ff = fn->fnct;
  while ( ff != 0 )
  {
    for ( i = 0; i < fn->outvars; i++ )
      if ( ff->ov[i] == '1' ) break;
    ff->valid = ( i == fn->outvars ) ? 0 : 1;
    ff = ff->next;
  }
  fn->DeleteInvalidTerms();
  gauge.CloseGauge();
}
//---------------------------------------------------------------------------
TBoolFnct * BlockCover( TBoolFnct * fn, Cov_Opts * covopts, TGauge gauge, int blocks )
{
  int i, j;
  CovMatrix * cm;
  TBoolFnct * outfn;
  TPfnct * ff1;
  CovLink * cl;
  CovSol * best;

  CleanInpAssignment( fn );

  outfn = fn->parent->CreateChildFnct();
  outfn->type = 0;

  for ( i = 0; i < blocks; i++ )
  {
    cm = CreateBlockMatrix( fn, i );
    best = ProcessCover( cm, covopts, gauge );
    cl = best->sol;
    while ( cl != 0 )
    {
      ff1 = outfn->AddTermToEnd( cl->term->iv, cl->term->ov, cl->term->dim );
      ff1->cov = DuplCov( cl->term->cov );
      for ( j = 0; j < fn->outvars; j++ )
        if ( fn->decomp->outs[j] != i ) ff1->ov[j] = '0';
      cl = cl->next;
    }
    DelLinks( best->sol );
    delete best;
    DeleteMatrix( cm );
  }
  return outfn;
}
//---------------------------------------------------------------------------
