SIMPLEX PRIMAL

/*******************************/
/*El Antri Abdellah            */
/* UDL                         */
/* Departement de mathématique */
/* sous licence GNU            */
/*******************************/




/* Comiler et executer sous unix system V */

#include<values.h>

typedef struct { int l,c; float tab[10][10];char nom[15];} matrice;
char car;
matrice A,b, CB, CN, XBB, XNB, N, B, I, NB, J;

int probleme_max = 1;

float val_abs(float value)
{
return ( value>0 ? value : -value);
}

void  produit (matrice m12, matrice m2, matrice *m3);
void  transpose(matrice m1, matrice *m2);
void  inverse(matrice m1, matrice *m2);
void  comatrice(matrice m1, matrice *m2);
void  lire(matrice *m);
void  afficher(matrice m);
float determinant(matrice m);
void  extraire (matrice m1, int l, int c, matrice *m2);
void  entrer_donnee();
void  initialisation();
void  extraire_col(matrice m1, int c, matrice *m2);
void  simplex();

/******************* Programmes principales ********************/

int main()
{
int i,j;
clrscr();

printf("\n Specifie le type du problème:");
scanf("%d",&probleme_max );
printf("Remplissage des matrices :\n");
initialisation();
simplex();
scanf("%c",car);
return 1;
}

/******************** fonctions et procedures *********************/
void initialisation()
{
int  i, j;

printf("\n Remplir la matrice A:");
strcpy(A.nom,"A");
lire(&A);afficher(A);
printf("\n Remplir la matrice b:");
strcpy(b.nom,"b");
lire(&b);
printf("\n Remplir la matrice CN:");
strcpy(CN.nom,"CN");
lire(&CN);
printf("\n Remplir la matrice I:");
strcpy(I.nom,"I");
lire(&I);
printf("\n Remplir la matrice NB:");
strcpy(NB.nom,"NB");
lire(&NB);

/* Calcule de CB */
strcpy(CB.nom,"CB");
for(i=0; i < I.c; i++)
       CB.tab[1][i] = 0;
CB.l = 1; CB.c = I.c;
afficher(CB);
scanf("%c",car);

/* Caclcul de la matrice N */
strcpy(N.nom,"N");
for(i=0; i < A.l; i++)
   for(j=0; j < NB.c; j++)
       N.tab[i][j] = A.tab[i][j];
N.l = A.l; N.c = NB.c;
afficher(N);
scanf("%c",car);
/* Caclcul de la matrice B */
strcpy(B.nom,"B");
for(i= 0; i < A.l ; i++)
   for(j=NB.c; j < A.c; j++)
       B.tab[i][j-NB.c] = A.tab[i][j];
B.l = A.l; B.c = I.c;
afficher(B);
scanf("%c",car);

/* Caclcul de la matrice XBB */
strcpy(XBB.nom,"XBB");
for(i=0; i < I.c; i++)
       XBB.tab[0][i] = b.tab[i][0];
XBB.l = 1; XBB.c = b.l;
afficher(XBB);
scanf("%c",car);

/* Caclcul de la matrice XNB */
strcpy(XNB.nom,"XNB");
for(i=0; i < NB.c; i++)
       XNB.tab[0][i] = 0;
XNB.l = 1; XNB.c = NB.c;
afficher(XNB);
scanf("%c",car);

/* Calcule de la matrice J*/
J.l = 1; J.c = NB.c;
for(i=0; i < NB.c; i++)
       J.tab[0][i] = NB.tab[0][i];
strcpy(J.nom,"J");
afficher(J);
scanf("%c",car);
}


void lire (matrice *m)
{
int i,j;
/* Lecture ligne par ligne */
printf("\nEntrez le nbr de ligne:");
scanf("%d",&m->l);
printf("\nEntrez le nbr de colonne:");
scanf("%d",&m->c);
for(i=0; i<m->l; i++)
    for(j=0; j<m->c; j++)
	{
	printf("\nEntrez l'element [%d][%d]:", i,j);
	scanf("%f",&m->tab[i][j]);
	}
}

void simplex()
{
int i,j, k;
float z[10], max;
float theta[10],t;
int entrant[10], i_e;
int sortant[10];
int fin;
matrice T1, T2, T3, IB;
matrice P;
int sol[10], i_s;

i_s = 0;
for(i=0; i < 10; i++) sol[i] = -1;
do
  {
  fin = 0; i_e = 0; max = 0; k=0;
  for(j = 0; j < I.c; j++)
      {
      theta[j] = MAXINT;
      }
  for(i = 0; i < NB.c; i++)
     {
     extraire_col(A,NB.tab[0][i]-1,&T1);
    inverse(B, &IB)   ;
     produit(CB,IB,&T2);
     produit(T2,T1,&T3);
     z[i] = T3.tab[0][0] - CN.tab[0][i];
     printf("\n z%d - c%d = %.2f",i+1,i+1,z[i]);
     scanf("%c",car);
     /* ProblŠme de max*/
    if(probleme_max)
      {
      if(z[i] < 0)
	{
	J.tab[0][i] = i;
	fin =1;
	}
	 else J.tab[0][i] = -1;
       }
       else /* Probleme de min */
	   {
	   if(z[i] > 0)
		{
		J.tab[0][i] = i;
		fin =1;
		}
		 else J.tab[0][i] = -1;
	   }
     }
  if(fin == 1)
    {
    for(j = 0; j < NB.c; j++)
      if(J.tab[0][j] >= 0)
	{
	extraire_col(A, J.tab[0][j], &T1);
	produit(IB,T1,&T2);
	sprintf(T2.nom,"inv(B)*a%d",j+1);
	afficher(T2);
	scanf("%c",car);
	i=0;
	while(i < B.l)
	{
	if(T2.tab[0][i] >= 0)break;
	i = i+1;
	}
	if(i == B.l)
	  {
	  printf("\n Pas de solution (Sol imborn‚e).");
	  scanf("%c",car);
	  exit(1);
	  }
	}

    for(i = 0; i < NB.c; i++)
       if(J.tab[0][i] >= 0 )
	 {
	 if(val_abs(z[i]) > max)
	   {
	   max = val_abs(z[i]);
	   }
	 }
    produit(IB,b, &T1);
    for(i = 0; i < NB.c; i++)
       {
       if(J.tab[0][i] >= 0 )
	 {
	 if(val_abs(z[i]) == max)
	   {
	   entrant[i_e] = i;
	   extraire_col(A, i, &T2);
	   produit(IB,T2,&P);
	   for(j = 0; j < I.c; j++)
	      if(P.tab[j][0] > 0)
		{
		t = T1.tab[j][0]/P.tab[j][0];
		printf("\n Theta%d = %.2f ",j+1,t);
		if(t < theta[i_e] )
		  {
		  theta[i_e]   = t ;
		  sortant[i_e] = j ;
		  }
		}
	   i_e = i_e + 1;
	   }
	 }
       }
    t=theta[0];
    for(i = 0; i<i_e; i++)
       if(theta[i] >= t)
	 {
	 k = i;
	 t = theta[i];
	 }
    printf("\n Sortant = X%.0f Entrant = X%.0f ",I.tab[0][sortant[k]], NB.tab[0][entrant[k]]);
    sol[entrant[k]] = sortant[k]; i_s = i_s + 1;
    scanf("%c",car);
    for(i = 0; i < B.l ; i++)
       {
       B.tab[i][sortant[k]] = N.tab[i][entrant[k]];
       }
    afficher(B)       ;scanf("%c",car);
    inverse(B, &IB)   ;afficher(IB);
    produit(IB,b,&XBB);afficher(XBB);
    i = I.tab[0][sortant[k]];
    I.tab[0][sortant[k]] = NB.tab[0][entrant[k]];
    NB.tab[0][entrant[k]] = i;
    i = CB.tab[0][sortant[k]];
    CB.tab[0][sortant[k]] = CN.tab[0][entrant[k]];
    CN.tab[0][entrant[k]] = i;
    afficher(I);afficher(NB); afficher(CB); afficher(CN);
    }else {
	  printf("\nL'Optimum est atteint.\n");
	  for(i = 0; i< NB.c; i++)
	     if(sol[i] == -1)printf("\n X%d = 0",i+1);
	       else printf("\n X%d = %.2f",i+1,XBB.tab[sol[i]][0]);
	  produit(CB,XBB,&T1);
	  printf(" \n ZB = %f", T1.tab[0][0]);
	  }

  } while(fin == 1);

}


void  extraire_col(matrice m1, int c, matrice *m2)
{
int i;

m2->l = m1.l; m2->c = 1;
for(i = 0; i < m1.l; i++)
   {
   m2->tab[i][0] = m1.tab[i][c];
   }
}

void afficher(matrice m)
{
int i,j;

printf("\n La matrice %s:",m.nom);
for(i=0; i<m.l; i++)
    {
    printf("\n|");
    for(j=0; j<m.c; j++)
	printf(" %.2f ",m.tab[i][j]);
    printf("|");
    }
}

void produit(matrice m1, matrice m2, matrice *m3)
{
int i,j,k;

if(m1.c != m2.l)
  {
  printf("\nErreur fatal: Produit de deux matrice incompatible");
  scanf("%c",car);
  exit(1);
  }

for(i=0; i<m2.l; i++)
  for(j=0; j<m2.c; j++)
     m3->tab[i][j] = 0;

m3->l = m1.l; m3->c = m2.c;

for(i=0; i<m1.l; i++)
  for(k=0; k<m2.l; k++)
     for(j=0; j<m1.c; j++)
	 m3->tab[i][k] = m3->tab[i][k] + m1.tab[i][j]*m2.tab[j][k];
}

void transpose(matrice m1, matrice *m2)
{
int i,j;

m2->l = m1.c;
m2->c = m1.l;
for(i=0; i<m1.l; i++)
    for(j=0; j<m1.c; j++)
	m2->tab[j][i] = m1.tab[i][j];
}

void extraire (matrice m1, int l, int c, matrice *m2)
{
int i,j;

for(i = l; i<m1.l; i++)
   for(j=0; j<m1.c; j++)
      m1.tab[i][j] = m1.tab[i+1][j];
for(j = c; j<m1.c; j++)
    for(i = 0; i<m1.l; i++)
	m1.tab[i][j] = m1.tab[i][j+1];
m2->l = m1.l -1;
m2->c = m1.c -1;
for(i=0; i<m2->l; i++)
   for(j=0; j<m2->l; j++)
       m2->tab[i][j] = m1.tab[i][j];
}

void comatrice(matrice m1, matrice *m2)
{
int i,j,p;
matrice m3;

m2->l = m1.l;
m2->c = m1.c;

for(i=0; i < m1.l; i++)
   for(j=0; j < m1.c; j++)
      {
      extraire(m1,i,j, &m3);
      if (((i+j) % 2) == 0) p = 1; else p = -1;
      m2->tab[i][j] = p*determinant(m3);
      }
}

float determinant(matrice m)
{
float d=0;
int k = 1;
int i,j;
matrice m2;

if( m.l != m.c)
    {
    printf("\n Erreur fatal: Determinant d'une matrice non carre.\n");
    scanf("%c",car);
    exit(1);
    }
    else
      {
      if( m.l == 1){ d =  m.tab[0][0];}
	else
	 for(j = 0; j<m.c; j++)
	   {
	   extraire(m,0,j,&m2);
	   d = d + k*m.tab[0][j]*determinant(m2);
	   k = -k;
	   }
      }
return d;
}

void inverse(matrice m1, matrice *m2)
{
float d;
int i,j;
matrice m,mm;

sprintf(m2->nom,"inv(%s)",m1.nom);
m2->l = m1.l; m2->c = m1.c;
d = determinant(m1);
comatrice(m1, &mm);
transpose(mm, &m);
for(i=0; i < m.l; i++)
   {
   for(j=0; j < m.c; j++)
      {
      m2->tab[i][j] = m.tab[i][j]/d;
      }
   }
}