simplex Class Reference

Linear Programming / Linear Optimization using Simplex - Algorithm. More...

#include <mpr_numeric.h>

Public Member Functions
	simplex (int rows, int cols)
	#rows should be >= m+2, #cols >= n+1 More...
	~simplex ()
BOOLEAN	mapFromMatrix (matrix m)
matrix	mapToMatrix (matrix m)
intvec *	posvToIV ()
intvec *	zrovToIV ()
void	compute ()

Data Fields
int	m
int	n
int	m1
int	m2
int	m3
int	icase
int *	izrov
int *	iposv
mprfloat **	LiPM

Private Member Functions
	simplex (const simplex &)
void	simp1 (mprfloat **a, int mm, int ll[], int nll, int iabf, int *kp, mprfloat *bmax)
void	simp2 (mprfloat **a, int n, int l2[], int nl2, int *ip, int kp, mprfloat *q1)
void	simp3 (mprfloat **a, int i1, int k1, int ip, int kp)

Private Attributes
int	LiPM_cols
int	LiPM_rows

Detailed Description

Linear Programming / Linear Optimization using Simplex - Algorithm.

On output, the tableau LiPM is indexed by two arrays of integers. ipsov[j] contains, for j=1..m, the number i whose original variable is now represented by row j+1 of LiPM. (left-handed vars in solution) (first row is the one with the objective function) izrov[j] contains, for j=1..n, the number i whose original variable x_i is now a right-handed variable, rep. by column j+1 of LiPM. These vars are all zero in the solution. The meaning of n<i<n+m1+m2 is the same as above.

Definition at line 194 of file mpr_numeric.h.

Constructor & Destructor Documentation

simplex() [1/2]

simplex::simplex	(	int	rows,
		int	cols
	)

#rows should be >= m+2, #cols >= n+1

Definition at line 972 of file mpr_numeric.cc.

   : LiPM_cols(cols), LiPM_rows(rows)
{
  int i;
 
  LiPM_rows=LiPM_rows+3;
  LiPM_cols=LiPM_cols+2;
 
  LiPM = (mprfloat **)omAlloc( LiPM_rows * sizeof(mprfloat *) );  // LP matrix
  for( i= 0; i < LiPM_rows; i++ )
  {
    // Mem must be allocated aligned, also for type double!
    LiPM[i] = (mprfloat *)omAlloc0Aligned( LiPM_cols * sizeof(mprfloat) );
  }
 
  iposv = (int *)omAlloc0( 2*LiPM_rows*sizeof(int) );
  izrov = (int *)omAlloc0( 2*LiPM_rows*sizeof(int) );
 
  m=n=m1=m2=m3=icase=0;
 
#ifdef mprDEBUG_ALL
  Print("LiPM size: %d, %d\n",LiPM_rows,LiPM_cols);
#endif
}

~simplex()

simplex::~simplex

(

)

Definition at line 997 of file mpr_numeric.cc.

{
  // clean up
  int i;
  for( i= 0; i < LiPM_rows; i++ )
  {
    omFreeSize( (void *) LiPM[i], LiPM_cols * sizeof(mprfloat) );
  }
  omFreeSize( (void *) LiPM, LiPM_rows * sizeof(mprfloat *) );
 
  omFreeSize( (void *) iposv, 2*LiPM_rows*sizeof(int) );
  omFreeSize( (void *) izrov, 2*LiPM_rows*sizeof(int) );
}

simplex() [2/2]

simplex::simplex ( const simplex &  )

private

Member Function Documentation

compute()

void simplex::compute

(

)

Definition at line 1095 of file mpr_numeric.cc.

{
  int i,ip,ir,is,k,kh,kp,m12,nl1,nl2;
  int *l1,*l2,*l3;
  mprfloat q1, bmax;
 
  if ( m != (m1+m2+m3) )
  {
    // error: bad input
    error(WarnS("simplex::compute: Bad input constraint counts!");)
    icase=-2;
    return;
  }
 
  l1= (int *) omAlloc0( (n+1) * sizeof(int) );
  l2= (int *) omAlloc0( (m+1) * sizeof(int) );
  l3= (int *) omAlloc0( (m+1) * sizeof(int) );
 
  nl1= n;
  for ( k=1; k<=n; k++ ) l1[k]=izrov[k]=k;
  nl2=m;
  for ( i=1; i<=m; i++ )
  {
    if ( LiPM[i+1][1] < 0.0 )
    {
      // error: bad input
      error(WarnS("simplex::compute: Bad input tableau!");)
      error(Warn("simplex::compute: in input Matrix row %d, column 1, value %f",i+1,LiPM[i+1][1]);)
      icase=-2;
      // free mem l1,l2,l3;
      omFreeSize( (void *) l3, (m+1) * sizeof(int) );
      omFreeSize( (void *) l2, (m+1) * sizeof(int) );
      omFreeSize( (void *) l1, (n+1) * sizeof(int) );
      return;
    }
    l2[i]= i;
    iposv[i]= n+i;
  }
  for ( i=1; i<=m2; i++) l3[i]= 1;
  ir= 0;
  if (m2+m3)
  {
    ir=1;
    for ( k=1; k <= (n+1); k++ )
    {
      q1=0.0;
      for ( i=m1+1; i <= m; i++ ) q1+= LiPM[i+1][k];
      LiPM[m+2][k]= -q1;
    }
 
    do
    {
      simp1(LiPM,m+1,l1,nl1,0,&kp,&bmax);
      if ( bmax <= SIMPLEX_EPS && LiPM[m+2][1] < -SIMPLEX_EPS )
      {
        icase= -1; // no solution found
        // free mem l1,l2,l3;
        omFreeSize( (void *) l3, (m+1) * sizeof(int) );
        omFreeSize( (void *) l2, (m+1) * sizeof(int) );
        omFreeSize( (void *) l1, (n+1) * sizeof(int) );
        return;
      }
      else if ( bmax <= SIMPLEX_EPS && LiPM[m+2][1] <= SIMPLEX_EPS )
      {
        m12= m1+m2+1;
        if ( m12 <= m )
        {
          for ( ip= m12; ip <= m; ip++ )
          {
            if ( iposv[ip] == (ip+n) )
            {
              simp1(LiPM,ip,l1,nl1,1,&kp,&bmax);
              if ( fabs(bmax) >= SIMPLEX_EPS)
                goto one;
            }
          }
        }
        ir= 0;
        --m12;
        if ( m1+1 <= m12 )
          for ( i=m1+1; i <= m12; i++ )
            if ( l3[i-m1] == 1 )
              for ( k=1; k <= n+1; k++ )
                LiPM[i+1][k] = -(LiPM[i+1][k]);
        break;
      }
      //#if DEBUG
      //print_bmat( a, m+2, n+3);
      //#endif
      simp2(LiPM,n,l2,nl2,&ip,kp,&q1);
      if ( ip == 0 )
      {
        icase = -1; // no solution found
        // free mem l1,l2,l3;
        omFreeSize( (void *) l3, (m+1) * sizeof(int) );
        omFreeSize( (void *) l2, (m+1) * sizeof(int) );
        omFreeSize( (void *) l1, (n+1) * sizeof(int) );
        return;
      }
    one: simp3(LiPM,m+1,n,ip,kp);
      // #if DEBUG
      // print_bmat(a,m+2,n+3);
      // #endif
      if ( iposv[ip] >= (n+m1+m2+1))
      {
        for ( k= 1; k <= nl1; k++ )
          if ( l1[k] == kp ) break;
        --nl1;
        for ( is=k; is <= nl1; is++ ) l1[is]= l1[is+1];
        ++(LiPM[m+2][kp+1]);
        for ( i= 1; i <= m+2; i++ ) LiPM[i][kp+1] = -(LiPM[i][kp+1]);
      }
      else
      {
        if ( iposv[ip] >= (n+m1+1) )
        {
          kh= iposv[ip]-m1-n;
          if ( l3[kh] )
          {
            l3[kh]= 0;
            ++(LiPM[m+2][kp+1]);
            for ( i=1; i<= m+2; i++ )
              LiPM[i][kp+1] = -(LiPM[i][kp+1]);
          }
        }
      }
      is= izrov[kp];
      izrov[kp]= iposv[ip];
      iposv[ip]= is;
    } while (ir);
  }
  /* end of phase 1, have feasible sol, now optimize it */
  loop
  {
    // #if DEBUG
    // print_bmat( a, m+1, n+5);
    // #endif
    simp1(LiPM,0,l1,nl1,0,&kp,&bmax);
    if (bmax <= /*SIMPLEX_EPS*/0.0)
    {
      icase=0; // finite solution found
      // free mem l1,l2,l3
      omFreeSize( (void *) l3, (m+1) * sizeof(int) );
      omFreeSize( (void *) l2, (m+1) * sizeof(int) );
      omFreeSize( (void *) l1, (n+1) * sizeof(int) );
      return;
    }
    simp2(LiPM,n,l2,nl2,&ip,kp,&q1);
    if (ip == 0)
    {
      //printf("Unbounded:");
      // #if DEBUG
      //       print_bmat( a, m+1, n+1);
      // #endif
      icase=1;                /* unbounded */
      // free mem
      omFreeSize( (void *) l3, (m+1) * sizeof(int) );
      omFreeSize( (void *) l2, (m+1) * sizeof(int) );
      omFreeSize( (void *) l1, (n+1) * sizeof(int) );
      return;
    }
    simp3(LiPM,m,n,ip,kp);
    is= izrov[kp];
    izrov[kp]= iposv[ip];
    iposv[ip]= is;
  }/*for ;;*/
}

mapFromMatrix()

BOOLEAN simplex::mapFromMatrix

(

matrix

m

)

Definition at line 1011 of file mpr_numeric.cc.

{
  int i,j;
//    if ( MATROWS( m ) > LiPM_rows ||  MATCOLS( m ) > LiPM_cols ) {
//      WarnS("");
//      return FALSE;
//    }
 
  number coef;
  for ( i= 1; i <= MATROWS( mm ); i++ )
  {
     for ( j= 1; j <= MATCOLS( mm ); j++ )
     {
        if ( MATELEM(mm,i,j) != NULL )
        {
           coef= pGetCoeff( MATELEM(mm,i,j) );
           if ( coef != NULL && !nIsZero(coef) )
              LiPM[i][j]= (double)(*(gmp_float*)coef);
           //#ifdef mpr_DEBUG_PROT
           //Print("%f ",LiPM[i][j]);
           //#endif
        }
     }
     //     PrintLn();
  }
 
  return TRUE;
}

mapToMatrix()

matrix simplex::mapToMatrix

(

matrix

m

)

Definition at line 1040 of file mpr_numeric.cc.

{
  int i,j;
//    if ( MATROWS( mm ) < LiPM_rows-3 ||  MATCOLS( m ) < LiPM_cols-2 ) {
//      WarnS("");
//      return NULL;
//    }
 
//Print(" %d x %d\n",MATROWS( mm ),MATCOLS( mm ));
 
  number coef;
  gmp_float * bla;
  for ( i= 1; i <= MATROWS( mm ); i++ )
  {
    for ( j= 1; j <= MATCOLS( mm ); j++ )
    {
       pDelete( &(MATELEM(mm,i,j)) );
       MATELEM(mm,i,j)= NULL;
//Print(" %3.0f ",LiPM[i][j]);
       if ( LiPM[i][j] != 0.0 )
       {
          bla= new gmp_float(LiPM[i][j]);
          coef= (number)bla;
          MATELEM(mm,i,j)= pOne();
          pSetCoeff( MATELEM(mm,i,j), coef );
       }
    }
//PrintLn();
  }
 
  return mm;
}

posvToIV()

intvec * simplex::posvToIV

(

)

Definition at line 1073 of file mpr_numeric.cc.

{
   int i;
   intvec * iv = new intvec( m );
   for ( i= 1; i <= m; i++ )
   {
      IMATELEM(*iv,i,1)= iposv[i];
   }
   return iv;
}

simp1()

void simplex::simp1 ( mprfloat **  a, int  mm, int  ll[], int  nll, int  iabf, int *  kp, mprfloat *  bmax  )

private

Definition at line 1263 of file mpr_numeric.cc.

{
  int k;
  mprfloat test;
 
  if( nll <= 0)
  {                        /* init'tion: fixed */
    *bmax = 0.0;
    return;
  }
  *kp=ll[1];
  *bmax=a[mm+1][*kp+1];
  for (k=2;k<=nll;k++)
  {
    if (iabf == 0)
    {
      test=a[mm+1][ll[k]+1]-(*bmax);
      if (test > 0.0)
      {
        *bmax=a[mm+1][ll[k]+1];
        *kp=ll[k];
      }
    }
    else
    {                        /* abs values: have fixed it */
      test=fabs(a[mm+1][ll[k]+1])-fabs(*bmax);
      if (test > 0.0)
      {
        *bmax=a[mm+1][ll[k]+1];
        *kp=ll[k];
      }
    }
  }
}

simp2()

void simplex::simp2 ( mprfloat **  a, int  n, int  l2[], int  nl2, int *  ip, int  kp, mprfloat *  q1  )

private

Definition at line 1298 of file mpr_numeric.cc.

{
  int k,ii,i;
  mprfloat qp,q0,q;
 
  *ip= 0;
  for ( i=1; i <= nl2; i++ )
  {
    if ( a[l2[i]+1][kp+1] < -SIMPLEX_EPS )
    {
      *q1= -a[l2[i]+1][1] / a[l2[i]+1][kp+1];
      *ip= l2[i];
      for ( i= i+1; i <= nl2; i++ )
      {
        ii= l2[i];
        if (a[ii+1][kp+1] < -SIMPLEX_EPS)
        {
          q= -a[ii+1][1] / a[ii+1][kp+1];
          if (q - *q1 < -SIMPLEX_EPS)
          {
            *ip=ii;
            *q1=q;
          }
          else if (q - *q1 < SIMPLEX_EPS)
          {
            for ( k=1; k<= nn; k++ )
            {
              qp= -a[*ip+1][k+1]/a[*ip+1][kp+1];
              q0= -a[ii+1][k+1]/a[ii+1][kp+1];
              if ( q0 != qp ) break;
            }
            if ( q0 < qp ) *ip= ii;
          }
        }
      }
    }
  }
}

simp3()

void simplex::simp3 ( mprfloat **  a, int  i1, int  k1, int  ip, int  kp  )

private

Definition at line 1337 of file mpr_numeric.cc.

{
  int kk,ii;
  mprfloat piv;
 
  piv= 1.0 / a[ip+1][kp+1];
  for ( ii=1; ii <= i1+1; ii++ )
  {
    if ( ii -1 != ip )
    {
      a[ii][kp+1] *= piv;
      for ( kk=1; kk <= k1+1; kk++ )
        if ( kk-1 != kp )
          a[ii][kk] -= a[ip+1][kk] * a[ii][kp+1];
    }
  }
  for ( kk=1; kk<= k1+1; kk++ )
    if ( kk-1 != kp ) a[ip+1][kk] *= -piv;
  a[ip+1][kp+1]= piv;
}

zrovToIV()

intvec * simplex::zrovToIV

(

)

Definition at line 1084 of file mpr_numeric.cc.

{
   int i;
   intvec * iv = new intvec( n );
   for ( i= 1; i <= n; i++ )
   {
      IMATELEM(*iv,i,1)= izrov[i];
   }
   return iv;
}

Field Documentation

icase

int simplex::icase

Definition at line 201 of file mpr_numeric.h.

iposv

int * simplex::iposv

Definition at line 203 of file mpr_numeric.h.

izrov

int* simplex::izrov

Definition at line 203 of file mpr_numeric.h.

LiPM

mprfloat** simplex::LiPM

Definition at line 205 of file mpr_numeric.h.

LiPM_cols

int simplex::LiPM_cols

private

Definition at line 225 of file mpr_numeric.h.

LiPM_rows

int simplex::LiPM_rows

private

Definition at line 225 of file mpr_numeric.h.

m

int simplex::m

Definition at line 198 of file mpr_numeric.h.

m1

int simplex::m1

Definition at line 200 of file mpr_numeric.h.

m2

int simplex::m2

Definition at line 200 of file mpr_numeric.h.

m3

int simplex::m3

Definition at line 200 of file mpr_numeric.h.

n

int simplex::n

Definition at line 199 of file mpr_numeric.h.

---

The documentation for this class was generated from the following files:

• kernel/numeric/mpr_numeric.h
• kernel/numeric/mpr_numeric.cc