sc_simplex_feasibility_fixprec.c File Reference

#include <stdlib.h>
#include <stdio.h>
#include <stdint.h>
#include <string.h>
#include <limits.h>
#include "boolean.h"
#include "arithmetique.h"
#include "linear_assert.h"
#include "vecteur.h"
#include "contrainte.h"
#include "sc.h"

Include dependency graph for sc_simplex_feasibility_fixprec.c:

Go to the source code of this file.

Data Structures
struct	hashtable_t
	this is already too much... More...

Macros
#define	DEBUG(code)   { }
	debug macros may be trigered with -DDEBUG{,1,2} More...
#define	DEBUG1(code)   { }
#define	DEBUG2(code)   { }
#define	DEBUG3(code)   { }
#define	DIMENSION   sc->dimension
#define	NUMERO   hashtable[h].numero
#define	SOLUBLE(N)   soluble=N;goto FINSIMPLEX ;
#define	CREVARVISIBLE   variables[compteur-3]=compteur-2;
#define	CREVARCACHEE
#define	tag(x)   /**printf(x); */
	for tracing macros after expansion: More...
#define	G(j, a, b)
	maybe most of them should be functions? More...
#define	GCD(j, a, b)
#define	GCDZZN(j, a, b)
#define	GCD_ZERO_CTRL(j, a, b)
#define	SIMPL(a, b)
	SIMPL normalizes rational a/b (b<>0): divides by gcd(a,b) and returns with b>0 note that there should be no arithmetic exceptions within this macro: (well, only uminus may have trouble for VALUE_MIN...) ??? a==0 ? then a = a/b = 0 and b = b/b = 1. More...
#define	SIMPLIFIE(f)
	SIMPLIFIE normalizes fraction f. More...
#define	AFF(x, y)   {x.num=y.num; x.den=y.den;} /**x=y should be ok:-) */
#define	AFF_PX(x, y)   {x->num=y.num; x->den=y.den;}
#define	INV(x, y)   {x.num=y.den; x.den=y.num;} /**x=1/y */
#define	INV_ZERO_CTRL(x, y)   {if (value_zero_p(y.num)) {fprintf(stderr,"ERROR : inverse of fraction zero !"); x.num = VALUE_ZERO; x.den = VALUE_ONE;} else {INV(x,y)}}
	??? value_zero_p(y.num)? Then x.num != VALUE_ZERO and x.den = VALUE_ZERO, it's not good at all. More...
#define	METINFINI(f)   {f.num=VALUE_MAX; f.den=VALUE_ONE;}
#define	MET_ZERO(f)   {f.num=VALUE_ZERO; f.den=VALUE_ONE;}
#define	MET_UN(f)   {f.num=VALUE_ONE; f.den=VALUE_ONE;}
#define	EGAL1(x)   (value_eq(x.num,x.den))
#define	EGAL0(x)   (value_zero_p(x.num))
#define	NUL(x)   (value_zero_p(x.num))
#define	NEGATIF(x)
#define	POSITIF(x)
#define	SUP1(x)
#define	EGAL_MACRO(x, y, mult)
#define	INF_MACRO(x, y, mult)   ((value_pos_p(mult(x.den,y.den)) && value_lt(mult(x.num,y.den),mult(x.den,y.num))) || (value_neg_p(mult(x.den,y.den)) && value_gt(mult(x.num,y.den),mult(x.den,y.num))))
	define INF_MACRO(x,y,mult) (value_lt(mult(x.num,y.den),mult(x.den,y.num))) DN: c'est pas assez pour la comparaison entre deux fractions, qu'est-ce qui se passe s'il y a seulement un denominateur negatif??? Ca donnera un resultat faux. More...
#define	DIV_MACRO(x, y, z, mult)
	computes x = simplify(y/z) More...
#define	MUL_MACRO(x, y, z, mult)
	computes x = simplify(y*z) More...
#define	SUB_MACRO(X, A, B, mult)
	computes X = simplify(A-B) More...
#define	FULL_PIVOT_MACRO_SIOUX(X, A, B, C, D, mult)
	computes X = A - B*C/D, trying to avoid arithmetic exceptions... More...
#define	FULL_PIVOT_MACRO_DIRECT(X, A, B, C, D, mult)
	computes X = A - B*C/D, but does not try to avoid arithmetic exceptions More...
#define	direct_p(v)   (value_lt(v,MAXVAL))
#define	FULL_PIVOT_MACRO(X, A, B, C, D, mult)
	computes X = A - B*C/D, with a switch to use SIOUX or DIRECT thae rationale for the actual condition is quite fuzzy. More...
#define	PARTIAL_PIVOT_MACRO_SIOUX(X, B, C, D, mult)
	idem if A==0 More...
#define	PARTIAL_PIVOT_MACRO_DIRECT(X, B, C, D, mult)
#define	PARTIAL_PIVOT_MACRO(X, B, C, D, mult)
#define	PIVOT_MACRO(X, A, B, C, D, mult)
	Pivot : x = a - b c / d mult is used for multiplying values. More...
#define	value_protected_mult(v, w)    value_protected_multiply(v,w,THROW(simplex_arithmetic_error))
	multiplies two Values of no arithmetic overflow, or throw exception. More...
#define	MULT(RES, A, B)   RES=value_mult(A,B)
	Version with and without arithmetic exceptions... More...
#define	MULTOFL(RES, A, B)   RES=value_protected_mult(A,B)
#define	DIV(x, y, z)   DIV_MACRO(x,y,z,value_mult)
#define	DIVOFL(x, y, z)   DIV_MACRO(x,y,z,value_protected_mult)
#define	MUL(x, y, z)   MUL_MACRO(x,y,z,value_mult)
#define	MULOFL(x, y, z)   MUL_MACRO(x,y,z,value_protected_mult)
#define	SUB(X, A, B)   SUB_MACRO(X,A,B,value_mult)
#define	SUBOFL(X, A, B)   SUB_MACRO(X,A,B,value_protected_mult)
#define	PIVOT(X, A, B, C, D)   PIVOT_MACRO(X,A,B,C,D,value_mult)
#define	PIVOTOFL(X, A, B, C, D)   PIVOT_MACRO(X,A,B,C,D,value_protected_mult)
#define	EGAL(x, y)   EGAL_MACRO(x,y,value_mult)
#define	EGALOFL(x, y)   EGAL_MACRO(x,y,value_protected_mult)
#define	INF(x, y)   INF_MACRO(x,y,value_mult)
#define	INFOFL(x, y)   INF_MACRO(x,y,value_protected_mult)

Enumerations
enum	{ PTR_NIL = INTPTR_MIN+767 , INFINI = INTPTR_MAX-767 , MAX_VAR = 1971 , MAXVAL = 24 }
	Hmmm. More...

Functions
void	frac_init (frac *f, int n)
void	frac_simpl (Value *a, Value *b)
void	frac_simplifie (frac *f)
	simplifie normalizes fraction f More...
void	frac_div (frac *x, frac y, frac z, bool ofl_ctrl)
void	frac_mul (frac *x, frac y, frac z, bool ofl_ctrl)
	computes x = simplify(y*z) More...
void	frac_sub (frac *X, frac A, frac B, bool ofl_ctrl)
void	full_pivot_sioux (frac *X, frac A, frac B, frac C, frac D, bool ofl_ctrl)
void	full_pivot_direct (frac *X, frac A, frac B, frac C, frac D, bool ofl_ctrl)
	computes X = A - B*C/D, but does not try to avoid arithmetic exceptions More...
void	full_pivot (frac *X, frac A, frac B, frac C, frac D, bool ofl_ctrl)
void	partial_pivot_sioux (frac *X, frac B, frac C, frac D, bool ofl_ctrl)
	idem if A==0 More...
void	partial_pivot_direct (frac *X, frac B, frac C, frac D, bool ofl_ctrl)
void	partial_pivot (frac *X, frac B, frac C, frac D, bool ofl_ctrl)
void	pivot (frac *X, frac A, frac B, frac C, frac D, bool ofl_ctrl)
static void	__attribute__ ((unused))
	For debugging: More...
static void	printfrac (frac x)
static int	hash (Variable s)
	dump_tableau More...
bool	sc_simplex_feasibility_ofl_ctrl_fixprec (Psysteme sc, int ofl_ctrl)
	fonction de calcul de la faisabilite' d'un systeme d'equations et d'inequations Auteur : Robert Mahl, Date : janvier 1994 More...

Variables
static int	NB_EQ = 0
	test du simplex : Si on compile grace a` "make simp" dans le repertoire /projects/C3/Linear/Development/polyedre/Tests alors on peut tester l'execution dans le meme directory en faisant : make test_simp More...
static int	NB_INEQ = 0
static frac	frac0 ={(Value)0,(Value)1,0}
static int	nbvariables
	Le nombre de variables visibles est : compteur-2 La i-eme variable visible a le numero : variables[i+1]=i (0 <= i < compteur-2) Le nombre de variables cachees est : nbvarables La i-eme variable cachee a le numero : variablescachees[i+1]=MAX_VAR+i-1 (0 <= i < nbvariables) More...
static int	variablescachees [MAX_VAR]
static int	variables [MAX_VAR]

Macro Definition Documentation

AFF

#define AFF	(		x,
			y
	)		{x.num=y.num; x.den=y.den;} /**x=y should be ok:-) */

Definition at line 224 of file sc_simplex_feasibility_fixprec.c.

AFF_PX

#define AFF_PX	(		x,
			y
	)		{x->num=y.num; x->den=y.den;}

Definition at line 225 of file sc_simplex_feasibility_fixprec.c.

CREVARCACHEE

#define CREVARCACHEE

Value:

                         { variablescachees[nbvariables]=nbvariables + MAX_VAR ; \
                         if (nbvariables ++ >= MAX_VAR) abort(); }

Definition at line 122 of file sc_simplex_feasibility_fixprec.c.

CREVARVISIBLE

#define CREVARVISIBLE   variables[compteur-3]=compteur-2;

Definition at line 121 of file sc_simplex_feasibility_fixprec.c.

DEBUG

#define DEBUG

(

code

)

{ }

debug macros may be trigered with -DDEBUG{,1,2}

Definition at line 61 of file sc_simplex_feasibility_fixprec.c.

DEBUG1

#define DEBUG1

(

code

)

{ }

Definition at line 68 of file sc_simplex_feasibility_fixprec.c.

DEBUG2

#define DEBUG2

(

code

)

{ }

Definition at line 75 of file sc_simplex_feasibility_fixprec.c.

DEBUG3

#define DEBUG3

(

code

)

{ }

Definition at line 82 of file sc_simplex_feasibility_fixprec.c.

DIMENSION

#define DIMENSION   sc->dimension

Definition at line 118 of file sc_simplex_feasibility_fixprec.c.

direct_p

#define direct_p

(

v

)

(value_lt(v,MAXVAL))

Definition at line 393 of file sc_simplex_feasibility_fixprec.c.

DIV

#define DIV	(		x,
			y,
			z
	)		DIV_MACRO(x,y,z,value_mult)

Definition at line 491 of file sc_simplex_feasibility_fixprec.c.

DIV_MACRO

#define DIV_MACRO	(		x,
			y,
			z,
			mult
	)

Value:

{tag("DIV_MACRO")                               \
    if (value_zero_p(y.num))                    \
    {                                           \
        MET_ZERO(x);                            \
    }                                           \
    else                                        \
    {                                   \
           if (value_zero_p(z.num))             \
           {                                    \
              fprintf(stderr,"ATTENTION : divided by zero number!");    \
              MET_ZERO(x);                      \
           }                                    \
           else                                 \
{                                               \
              x.num=mult(y.num,z.den);          \
              x.den=mult(y.den,z.num);          \
              SIMPLIFIE(x);                     \
}                                               \
    }                                           \
}

computes x = simplify(y/z)

define DIV_MACRO(x,y,z,mult) \ {tag("DIV_MACRO") \ if (value_zero_p(y.num)) \ { \ MET_ZERO(x); \ } \ else \ { \ x.num=mult(y.num,z.den); \ x.den=mult(y.den,z.num); \ SIMPLIFIE(x); \ } \ } DN: This macro DIV_MACRO doesn't test if z = 0, then x = y/z means x.num = y.num*z.den and x.den = 0. We better add the test : if z = 0 then x.num = 0 and x.den = 1 We try to avoid the denominator equal to 0.

Definition at line 292 of file sc_simplex_feasibility_fixprec.c.

DIVOFL

#define DIVOFL	(		x,
			y,
			z
	)		DIV_MACRO(x,y,z,value_protected_mult)

Definition at line 492 of file sc_simplex_feasibility_fixprec.c.

EGAL

#define EGAL	(		x,
			y
	)		EGAL_MACRO(x,y,value_mult)

Definition at line 503 of file sc_simplex_feasibility_fixprec.c.

EGAL0

#define EGAL0

(

x

)

(value_zero_p(x.num))

Definition at line 240 of file sc_simplex_feasibility_fixprec.c.

EGAL1

#define EGAL1

(

x

)

(value_eq(x.num,x.den))

Definition at line 239 of file sc_simplex_feasibility_fixprec.c.

EGAL_MACRO

#define EGAL_MACRO	(		x,
			y,
			mult
	)

Value:

  ((value_zero_p(x.num) && value_zero_p(y.num)) ||              \
   (value_notzero_p(x.den) && value_notzero_p(y.den) &&         \
    value_eq(mult(x.num,y.den),mult(x.den,y.num))))

Definition at line 255 of file sc_simplex_feasibility_fixprec.c.

EGALOFL

#define EGALOFL	(		x,
			y
	)		EGAL_MACRO(x,y,value_protected_mult)

Definition at line 504 of file sc_simplex_feasibility_fixprec.c.

FULL_PIVOT_MACRO

#define FULL_PIVOT_MACRO	(		X,
			A,
			B,
			C,
			D,
			mult
	)

Value:

{ tag("FULL_PIVOT")                                                     \
    if (direct_p(A.den) && direct_p(B.den) &&                           \
        direct_p(C.den) && direct_p(value_abs(D.num)))                  \
    {                                                                   \
        FULL_PIVOT_MACRO_DIRECT(X,A,B,C,D,mult);                        \
    }                                                                   \
    else                                                                \
    {                                                                   \
        FULL_PIVOT_MACRO_SIOUX(X,A,B,C,D,mult);                         \
    }                                                                   \
}

computes X = A - B*C/D, with a switch to use SIOUX or DIRECT thae rationale for the actual condition is quite fuzzy.

Definition at line 398 of file sc_simplex_feasibility_fixprec.c.

FULL_PIVOT_MACRO_DIRECT

#define FULL_PIVOT_MACRO_DIRECT	(		X,
			A,
			B,
			C,
			D,
			mult
	)

Value:

{                                                                         \
    Value v; tag("FULL_PIVOT_DIRECT")                                     \
    X.num = mult(A.num,B.den);                                            \
    X.num = mult(X.num,C.den);                                            \
    X.num = mult(X.num,D.num);                                            \
    v = mult(A.den,B.num);                                                \
    v = mult(v,C.num);                                                    \
    v = mult(v,D.den);                                                    \
    value_substract(X.num,v);                                             \
    X.den = mult(A.den,B.den);                                            \
    X.den = mult(X.den,C.den);                                            \
    X.den = mult(X.den,D.num);                                            \
    SIMPLIFIE(X);                                                         \
}

computes X = A - B*C/D, but does not try to avoid arithmetic exceptions

Definition at line 377 of file sc_simplex_feasibility_fixprec.c.

FULL_PIVOT_MACRO_SIOUX

#define FULL_PIVOT_MACRO_SIOUX	(		X,
			A,
			B,
			C,
			D,
			mult
	)

Value:

{                                                                       \
    frac u,v,w; tag("FULL_PIVOT_SIOUX")                                 \
    AFF(u,B); AFF(v,C); INV_ZERO_CTRL(w,D); /**u*v*w == B*C/D */                        \
    SIMPL(u.num,v.den); SIMPL(u.num,w.den);                             \
    SIMPL(v.num,u.den); SIMPL(v.num,w.den);                             \
    SIMPL(w.num,u.den); SIMPL(w.num,v.den);                             \
    u.num = mult(u.num,v.num); /**u*=v */                               \
    u.den = mult(u.den,v.den);                                          \
    u.num = mult(u.num,w.num); /**u*=w */                               \
    u.den = mult(u.den,w.den);                                          \
    SUB_MACRO(X,A,u,mult);                                              \
}

computes X = A - B*C/D, trying to avoid arithmetic exceptions...

Definition at line 361 of file sc_simplex_feasibility_fixprec.c.

G

#define G	(		j,
			a,
			b
	)

Value:

{tag("G")                                               \
    j=b;                                                \
    if (value_gt(b,VALUE_ONE))                          \
    {                                                   \
        Value i=a, k;                                   \
        while(value_notzero_p(k=value_mod(i,j)))        \
            i=j, j=k;                                   \
    }                                                   \
    if (value_neg_p(j))                                 \
        value_oppose(j);                                \
}

maybe most of them should be functions?

G computes j=gcd(a,b) assuming b>0 and better with a>b. there can be no artihmetic errors.

Definition at line 136 of file sc_simplex_feasibility_fixprec.c.

GCD

#define GCD	(		j,
			a,
			b
	)

Value:

{tag("GCD")                                     \
    if (value_gt(a,b))                          \
      { G(j,a,b); }                             \
    else                                        \
      { G(j,b,a); }                             \
}

Definition at line 149 of file sc_simplex_feasibility_fixprec.c.

GCD_ZERO_CTRL

#define GCD_ZERO_CTRL	(		j,
			a,
			b
	)

Value:

{tag("GCD_ZERO_CTRL")                                   \
    if ((value_notzero_p(a)) && (value_notzero_p(b)))   \
        { GCDZZN(j,a,b);}                               \
    else \
        {fprintf(stderr,"********************************************************************Error : GCD of number zero !");    \
        j = (VALUE_ONE);}                       \
}

Definition at line 181 of file sc_simplex_feasibility_fixprec.c.

GCDZZN

#define GCDZZN	(		j,
			a,
			b
	)

Value:

{tag("GCDZZN")                                          \
    Value i,k;                                          \
    i = (value_absolute(a)), k = (value_absolute(b));   \
    while(value_notzero_p(j=value_mod(i,k)))    \
            i=k, k=j;                                   \
    j = k;                                              \
}

Definition at line 172 of file sc_simplex_feasibility_fixprec.c.

INF

#define INF	(		x,
			y
	)		INF_MACRO(x,y,value_mult)

Definition at line 506 of file sc_simplex_feasibility_fixprec.c.

INF_MACRO

#define INF_MACRO	(		x,
			y,
			mult
	)		((value_pos_p(mult(x.den,y.den)) && value_lt(mult(x.num,y.den),mult(x.den,y.num))) || (value_neg_p(mult(x.den,y.den)) && value_gt(mult(x.num,y.den),mult(x.den,y.num))))

define INF_MACRO(x,y,mult) (value_lt(mult(x.num,y.den),mult(x.den,y.num))) DN: c'est pas assez pour la comparaison entre deux fractions, qu'est-ce qui se passe s'il y a seulement un denominateur negatif??? Ca donnera un resultat faux.

a/b < c/d <=> if b*d > 0 then a*d < b*c else a*d < b*c

Definition at line 266 of file sc_simplex_feasibility_fixprec.c.

INFOFL

#define INFOFL	(		x,
			y
	)		INF_MACRO(x,y,value_protected_mult)

Definition at line 507 of file sc_simplex_feasibility_fixprec.c.

INV

#define INV	(		x,
			y
	)		{x.num=y.den; x.den=y.num;} /**x=1/y */

Definition at line 226 of file sc_simplex_feasibility_fixprec.c.

INV_ZERO_CTRL

#define INV_ZERO_CTRL	(		x,
			y
	)		{if (value_zero_p(y.num)) {fprintf(stderr,"ERROR : inverse of fraction zero !"); x.num = VALUE_ZERO; x.den = VALUE_ONE;} else {INV(x,y)}}

??? value_zero_p(y.num)? Then x.num != VALUE_ZERO and x.den = VALUE_ZERO, it's not good at all.

(assuming y.den != VALUE_ZERO) This means : test if y = 0 then x = 0 else x = 1/y Change in line 286 : DN.

Definition at line 232 of file sc_simplex_feasibility_fixprec.c.

MET_UN

#define MET_UN

(

f

)

{f.num=VALUE_ONE; f.den=VALUE_ONE;}

Definition at line 237 of file sc_simplex_feasibility_fixprec.c.

MET_ZERO

#define MET_ZERO

(

f

)

{f.num=VALUE_ZERO; f.den=VALUE_ONE;}

Definition at line 236 of file sc_simplex_feasibility_fixprec.c.

METINFINI

#define METINFINI

(

f

)

{f.num=VALUE_MAX; f.den=VALUE_ONE;}

Definition at line 235 of file sc_simplex_feasibility_fixprec.c.

MUL

#define MUL	(		x,
			y,
			z
	)		MUL_MACRO(x,y,z,value_mult)

Definition at line 494 of file sc_simplex_feasibility_fixprec.c.

MUL_MACRO

#define MUL_MACRO	(		x,
			y,
			z,
			mult
	)

Value:

{tag("MUL_MACRO")                                       \
    if(value_zero_p(y.num) || value_zero_p(z.num))      \
        MET_ZERO(x)                                     \
    else                                                \
    {                                                   \
        x.num=mult(y.num,z.num);                        \
        x.den=mult(y.den,z.den);                        \
        SIMPLIFIE(x);                                   \
    }                                                   \
}

computes x = simplify(y*z)

Definition at line 317 of file sc_simplex_feasibility_fixprec.c.

MULOFL

#define MULOFL	(		x,
			y,
			z
	)		MUL_MACRO(x,y,z,value_protected_mult)

Definition at line 495 of file sc_simplex_feasibility_fixprec.c.

MULT

#define MULT	(		RES,
			A,
			B
	)		RES=value_mult(A,B)

Version with and without arithmetic exceptions...

Definition at line 488 of file sc_simplex_feasibility_fixprec.c.

MULTOFL

#define MULTOFL	(		RES,
			A,
			B
	)		RES=value_protected_mult(A,B)

Definition at line 489 of file sc_simplex_feasibility_fixprec.c.

NEGATIF

#define NEGATIF

(

x

)

Value:

  ((value_neg_p(x.num) && value_pos_p(x.den)) ||        \
   (value_pos_p(x.num) && value_neg_p(x.den)))

Definition at line 243 of file sc_simplex_feasibility_fixprec.c.

NUL

#define NUL

(

x

)

(value_zero_p(x.num))

Definition at line 241 of file sc_simplex_feasibility_fixprec.c.

NUMERO

#define NUMERO   hashtable[h].numero

Definition at line 119 of file sc_simplex_feasibility_fixprec.c.

PARTIAL_PIVOT_MACRO

#define PARTIAL_PIVOT_MACRO	(		X,
			B,
			C,
			D,
			mult
	)

Value:

{                                                               \
    if (direct_p(B.den) && direct_p(C.den) && direct_p(D.num))  \
    {                                                           \
        PARTIAL_PIVOT_MACRO_DIRECT(X,B,C,D,mult);               \
    }                                                           \
    else                                                        \
    {                                                           \
        PARTIAL_PIVOT_MACRO_SIOUX(X,B,C,D,mult);                \
    }                                                           \
}

Definition at line 431 of file sc_simplex_feasibility_fixprec.c.

PARTIAL_PIVOT_MACRO_DIRECT

#define PARTIAL_PIVOT_MACRO_DIRECT	(		X,
			B,
			C,
			D,
			mult
	)

Value:

{ tag("PARTIAL_PIVOT_DIRECT")                           \
    X.num = mult(B.num,C.num);                          \
    X.num = mult(X.num,D.den);                          \
    value_oppose(X.num);                                \
    X.den = mult(B.den,C.den);                          \
    X.den = mult(X.den,D.num);                          \
    SIMPLIFIE(X);                                       \
}

Definition at line 421 of file sc_simplex_feasibility_fixprec.c.

PARTIAL_PIVOT_MACRO_SIOUX

#define PARTIAL_PIVOT_MACRO_SIOUX	(		X,
			B,
			C,
			D,
			mult
	)

Value:

{ tag("PARTIAL_PIVOT_SIOUX")                            \
    frac u;                                             \
    MUL_MACRO(u,B,C,mult); /**u=simplify(b*c) */        \
    DIV_MACRO(X,u,D,mult); /**x=simplify(u/d) */        \
    value_oppose(X.num);   /**x=-x */                   \
}

idem if A==0

Definition at line 413 of file sc_simplex_feasibility_fixprec.c.

PIVOT

#define PIVOT	(		X,
			A,
			B,
			C,
			D
	)		PIVOT_MACRO(X,A,B,C,D,value_mult)

Definition at line 500 of file sc_simplex_feasibility_fixprec.c.

PIVOT_MACRO

#define PIVOT_MACRO	(		X,
			A,
			B,
			C,
			D,
			mult
	)

Value:

{ if (value_zero_p(D.num)) fprintf(stderr,"division of zero!!!");             \
    DEBUG3(fprintf(stdout, "pivot on: ");                                     \
           printfrac(A); printfrac(B); printfrac(C); printfrac(D));           \
   if (value_zero_p(A.num))/**a==0? */                                        \
   {                                                                          \
       if (value_zero_p(B.num) || value_zero_p(C.num) || value_zero_p(D.den)) \
           { MET_ZERO(X);}                                                    \
       else /**b*c/d != 0, calculons! */                                      \
           { PARTIAL_PIVOT_MACRO(X,B,C,D,mult);}                              \
   }                                                                          \
   else /**a!=0 */                                                            \
      if (value_zero_p(B.num) || value_zero_p(C.num) || value_zero_p(D.den))  \
          { AFF(X,A);}                                                        \
      else /** b*c/d != 0, calculons! */                                      \
          if (value_one_p(D.num) && value_one_p(A.den) &&                     \
              value_one_p(B.den) && value_one_p(C.den))                       \
          { /**no den to compute */                                           \
              Value v = mult(B.num,C.num);                                    \
              v = mult(v,D.den);                                              \
              X.num=value_minus(A.num,v);                                     \
              X.den=VALUE_ONE;                                                \
          }                                                                   \
          else /**well, we must compute the full formula! */                  \
              { FULL_PIVOT_MACRO(X,A,B,C,D,mult);}                            \
    DEBUG3(fprintf(stdout, " = "); printfrac(X); fprintf(stdout, "\n"));      \
}

Pivot : x = a - b c / d mult is used for multiplying values.

the macro has changed a lot, for indentation and so... FC. DN: Why don't we test d.num = 0 ? (meanwhile, we do test d.den = 0)

Definition at line 449 of file sc_simplex_feasibility_fixprec.c.

PIVOTOFL

#define PIVOTOFL	(		X,
			A,
			B,
			C,
			D
	)		PIVOT_MACRO(X,A,B,C,D,value_protected_mult)

Definition at line 501 of file sc_simplex_feasibility_fixprec.c.

POSITIF

#define POSITIF

(

x

)

Value:

  ((value_pos_p(x.num) && value_pos_p(x.den)) ||        \
   (value_neg_p(x.num) && value_neg_p(x.den)))

Definition at line 247 of file sc_simplex_feasibility_fixprec.c.

SIMPL

#define SIMPL	(		a,
			b
	)

Value:

{                                                       \
    if (value_notone_p(b) && value_notone_p(a) &&       \
        value_notmone_p(b) && value_notmone_p(a))       \
    {                                                   \
        register Value i=a, j=b, k;                     \
        while (value_notzero_p(k=value_mod(i,j)))       \
            i=j, j=k;                                   \
        value_division(a,j), value_division(b,j);       \
    }                                                   \
    if (value_neg_p(b))                                 \
        value_oppose(a), value_oppose(b);               \
}

SIMPL normalizes rational a/b (b<>0): divides by gcd(a,b) and returns with b>0 note that there should be no arithmetic exceptions within this macro: (well, only uminus may have trouble for VALUE_MIN...) ??? a==0 ? then a = a/b = 0 and b = b/b = 1.

Definition at line 197 of file sc_simplex_feasibility_fixprec.c.

SIMPLIFIE

#define SIMPLIFIE

(

f

)

Value:

{tag("SIMPLIFIE")                               \
     if (value_zero_p(f.den))                    \
       { MET_ZERO(f);}                           \
     else {                                      \
       if (value_zero_p(f.num))                 \
         f.den = VALUE_ONE;                     \
     else					\
       SIMPL(f.num,f.den);      }               \
}

SIMPLIFIE normalizes fraction f.

Definition at line 213 of file sc_simplex_feasibility_fixprec.c.

SOLUBLE

#define SOLUBLE

(

N

)

soluble=N;goto FINSIMPLEX ;

Definition at line 120 of file sc_simplex_feasibility_fixprec.c.

SUB

#define SUB	(		X,
			A,
			B
	)		SUB_MACRO(X,A,B,value_mult)

Definition at line 497 of file sc_simplex_feasibility_fixprec.c.

SUB_MACRO

#define SUB_MACRO	(		X,
			A,
			B,
			mult
	)

Value:

{ tag("SUB_MACRO")                                                            \
    if (value_zero_p(A.num))                                                  \
        X.num = value_uminus(B.num),                                          \
        X.den = B.den;                                                        \
    else if (value_zero_p(B.num))                                             \
    { AFF(X, A); }                                                            \
    else if (value_eq(A.den,B.den))                                           \
    {                                                                         \
        X.num = value_minus(A.num,B.num);                                     \
        X.den = A.den;                                                        \
        if (value_notone_p(A.den))                                            \
            { SIMPLIFIE(X);}                                                  \
    }                                                                         \
    else /**must compute the stuff: */                                        \
    {                                                                         \
        Value ad=A.den, bd=B.den, gd, v;                                      \
        GCD_ZERO_CTRL(gd,ad,bd);                                                              \
        if (value_notone_p(gd)) value_division(ad,gd), value_division(bd,gd); \
        X.num = mult(A.num,bd);                                               \
        v = mult(B.num,ad);                                                   \
        value_substract(X.num,v);                                             \
        v = mult(ad,bd);                                                      \
        X.den = mult(v,gd);                                                   \
        SIMPLIFIE(X);                                                         \
    }                                                                         \
}

computes X = simplify(A-B)

Definition at line 331 of file sc_simplex_feasibility_fixprec.c.

SUBOFL

#define SUBOFL	(		X,
			A,
			B
	)		SUB_MACRO(X,A,B,value_protected_mult)

Definition at line 498 of file sc_simplex_feasibility_fixprec.c.

SUP1

#define SUP1

(

x

)

Value:

  ((value_pos_p(x.num) && value_pos_p(x.den) &&  value_gt(x.num,x.den)) || \
   (value_neg_p(x.num) && value_neg_p(x.den) &&  value_gt(x.den,x.num)))

Definition at line 251 of file sc_simplex_feasibility_fixprec.c.

tag

#define tag

(

x

)

/**printf(x); */

for tracing macros after expansion:

Definition at line 127 of file sc_simplex_feasibility_fixprec.c.

value_protected_mult

#define value_protected_mult	(		v,
			w
	)		value_protected_multiply(v,w,THROW(simplex_arithmetic_error))

multiplies two Values of no arithmetic overflow, or throw exception.

this version is local to the simplex. note that under some defined macros value_mult can expand to value_protected_mult, which would be ok.

Definition at line 483 of file sc_simplex_feasibility_fixprec.c.

Enumeration Type Documentation

anonymous enum

anonymous enum

Hmmm.

To be compatible with some weird old 16-bit constants... RK

Enumerator
PTR_NIL
INFINI
MAX_VAR
MAXVAL	nombre max de variables seuil au dela duquel on se mefie d'un overflow

Definition at line 104 of file sc_simplex_feasibility_fixprec.c.

     {
  PTR_NIL = INTPTR_MIN+767,
  INFINI = INTPTR_MAX-767,
  MAX_VAR = 1971, /* nombre max de variables */
  /* seuil au dela duquel on se mefie d'un overflow
   */
#if defined(LINEAR_VALUE_IS_LONGLONG)
  MAXVAL = 576,
#else
  MAXVAL = 24  
#endif
};

Function Documentation

__attribute__()

static void __attribute__ ( (unused)  )

static

For debugging:

Definition at line 812 of file sc_simplex_feasibility_fixprec.c.

                                        {
  int i;
  for(i=0;i<MAX_VAR;i++)
    if(VARIABLE_DEFINED_P(hashtable[i].nom))
      printf("%s %d ", (char *) hashtable[i].nom, hashtable[i].numero),
        print_Value(hashtable[i].val),
        printf("\n");
}

References MAX_VAR, print_Value(), printf(), and VARIABLE_DEFINED_P.

Here is the call graph for this function:

frac_div()

void frac_div	(	frac *	x,
		frac	y,
		frac	z,
		bool	ofl_ctrl
	)

Definition at line 567 of file sc_simplex_feasibility_fixprec.c.

{tag("FRAC_DIV")
    if (value_zero_p(y.num))
      {
        MET_ZERO((*x));
      }
    else
      {
        if (value_zero_p(z.num))
          {
            fprintf(stderr,"ATTENTION : divided by zero number!");
            MET_ZERO((*x));
          }
        else
          {
            if (ofl_ctrl == FWD_OFL_CTRL) {
            x->num=value_protected_mult(y.num,z.den);
            x->den=value_protected_mult(y.den,z.num);
            }
            else {
              x->num=value_mult(y.num,z.den);
              x->den=value_mult(y.den,z.num);
            }
            frac_simplifie(x);
          }
      }
}

References frac::den, fprintf(), frac_simplifie(), FWD_OFL_CTRL, MET_ZERO, frac::num, tag, value_mult, value_protected_mult, value_zero_p, and x.

Referenced by partial_pivot_sioux(), and sc_simplex_feasibility_ofl_ctrl_fixprec().

Here is the call graph for this function:

Here is the caller graph for this function:

frac_init()

void frac_init	(	frac *	f,
		int	n
	)

Definition at line 523 of file sc_simplex_feasibility_fixprec.c.

{
  int i;
  for (i=0;i<n;i++) {
    f[i].num =  (Value)0;
    f[i].den =  (Value)1;
    f[i].numero =  0;
  }
}

References f().

Referenced by sc_simplex_feasibility_ofl_ctrl_fixprec().

Here is the call graph for this function:

Here is the caller graph for this function:

frac_mul()

void frac_mul	(	frac *	x,
		frac	y,
		frac	z,
		bool	ofl_ctrl
	)

computes x = simplify(y*z)

Definition at line 598 of file sc_simplex_feasibility_fixprec.c.

{tag("FRAC_MUL")
    if(value_zero_p(y.num) || value_zero_p(z.num))
      MET_ZERO((*x))
      else
        {
          if (ofl_ctrl == FWD_OFL_CTRL) {
          x->num=value_protected_mult(y.num,z.num);
          x->den=value_protected_mult(y.den,z.den);
          }
          else {
            x->num=value_mult(y.num,z.num);
            x->den=value_mult(y.den,z.den);
          }
          frac_simplifie(x);
        }
}

References frac::den, frac_simplifie(), FWD_OFL_CTRL, MET_ZERO, frac::num, tag, value_mult, value_protected_mult, value_zero_p, and x.

Referenced by partial_pivot_sioux().

Here is the call graph for this function:

Here is the caller graph for this function:

frac_simpl()

void frac_simpl	(	Value *	a,
		Value *	b
	)

Definition at line 533 of file sc_simplex_feasibility_fixprec.c.

{
  if (value_notone_p(*b) && value_notone_p(*a) &&
      value_notmone_p(*b) && value_notmone_p(*a))
    {
      register long long int lli = ABS(VALUE_TO_LONG(*a)), llj= ABS(VALUE_TO_LONG(*b)),k;
      if (lli<llj)  {
        long long int tmp = lli;
        lli=llj;llj=tmp;}
 
      while ((k=lli%llj) !=0)
        lli=llj, llj=k;
      value_division(*a,llj), value_division(*b,llj);
    }
  if (value_neg_p(*b))
    value_oppose(*a), value_oppose(*b);
}

References ABS, value_division, value_neg_p, value_notmone_p, value_notone_p, value_oppose, and VALUE_TO_LONG.

Referenced by frac_simplifie(), and full_pivot_sioux().

Here is the caller graph for this function:

frac_simplifie()

void frac_simplifie

(

frac *

f

)

simplifie normalizes fraction f

Definition at line 553 of file sc_simplex_feasibility_fixprec.c.

{
  tag("frac_simplifie")
    if (value_zero_p(f->den))
      {
        MET_ZERO((*f));}
    else {
      if (value_zero_p(f->num))
        f->den = VALUE_ONE;
      else
        frac_simpl(&(f->num),&(f->den));
    }
}

References f(), frac_simpl(), MET_ZERO, tag, VALUE_ONE, and value_zero_p.

Referenced by frac_div(), frac_mul(), frac_sub(), full_pivot_direct(), partial_pivot_direct(), and sc_simplex_feasibility_ofl_ctrl_fixprec().

Here is the call graph for this function:

Here is the caller graph for this function:

frac_sub()

void frac_sub	(	frac *	X,
		frac	A,
		frac	B,
		bool	ofl_ctrl
	)

must compute the stuff:

Definition at line 615 of file sc_simplex_feasibility_fixprec.c.

{ tag("FRAC_SUB")
    if (value_zero_p(A.num))
      X->num = value_uminus(B.num),
        X->den = B.den;
    else if (value_zero_p(B.num))
      { AFF_PX(X, A); }
    else if (value_eq(A.den,B.den))
      {
        X->num = value_minus(A.num,B.num);
        X->den = A.den;
        if (value_notone_p(A.den))
          { frac_simplifie(X);}
      }
    else /* must compute the stuff: */
      {
        Value ad=A.den, bd=B.den, gd, v;
        GCD_ZERO_CTRL(gd,ad,bd);
        if (value_notone_p(gd)) value_division(ad,gd), value_division(bd,gd);
        if (ofl_ctrl == FWD_OFL_CTRL) {
          X->num = value_protected_mult(A.num,bd);
          v = value_protected_mult(B.num,ad);
          value_substract(X->num,v);
          v =value_protected_mult(ad,bd);
          X->den = value_protected_mult(v,gd);
        }
        else {
          X->num = value_mult(A.num,bd);
          v = value_mult(B.num,ad);
          value_substract(X->num,v);
          v =value_mult(ad,bd);
          X->den = value_mult(v,gd);
        }
        frac_simplifie(X);
      }
}

References AFF_PX, B, frac_simplifie(), FWD_OFL_CTRL, GCD_ZERO_CTRL, tag, value_division, value_eq, value_minus, value_mult, value_notone_p, value_protected_mult, value_substract, value_uminus, value_zero_p, and X.

Referenced by full_pivot_sioux().

Here is the call graph for this function:

Here is the caller graph for this function:

full_pivot()

void full_pivot	(	frac *	X,
		frac	A,
		frac	B,
		frac	C,
		frac	D,
		bool	ofl_ctrl
	)

Definition at line 710 of file sc_simplex_feasibility_fixprec.c.

{ tag("FULL_PIVOT")
 
    if (direct_p(A.den) && direct_p(B.den) &&
        direct_p(C.den) && direct_p(value_abs(D.num)))
      {
        full_pivot_direct(X,A,B,C,D,ofl_ctrl);
      }
    else
      {
        full_pivot_sioux(X,A,B,C,D,ofl_ctrl);
      }
}

References B, D, frac::den, direct_p, full_pivot_direct(), full_pivot_sioux(), tag, value_abs, and X.

Referenced by pivot().

Here is the call graph for this function:

Here is the caller graph for this function:

full_pivot_direct()

void full_pivot_direct	(	frac *	X,
		frac	A,
		frac	B,
		frac	C,
		frac	D,
		bool	ofl_ctrl
	)

computes X = A - B*C/D, but does not try to avoid arithmetic exceptions

Definition at line 681 of file sc_simplex_feasibility_fixprec.c.

{
  Value v; tag("FULL_PIVOT_DIRECT") 
             if (ofl_ctrl == FWD_OFL_CTRL) {
               X->num = value_protected_mult(A.num,B.den);
               X->num = value_protected_mult(X->num,C.den);
               X->num = value_protected_mult(X->num,D.num);
               v = value_protected_mult(A.den,B.num);
               v = value_protected_mult(v,C.num);
               v = value_protected_mult(v,D.den);
               value_substract(X->num,v);
               X->den = value_protected_mult(A.den,B.den);
               X->den = value_protected_mult(X->den,C.den);
               X->den = value_protected_mult(X->den,D.num);
             }
             else {
               X->num = value_mult(A.num,B.den);
               X->num = value_mult(X->num,C.den);
               X->num = value_mult(X->num,D.num);
               v = value_mult(A.den,B.num);
               v = value_mult(v,C.num);
               v = value_mult(v,D.den);
               value_substract(X->num,v);
               X->den = value_mult(A.den,B.den);
               X->den = value_mult(X->den,C.den);
               X->den = value_mult(X->den,D.num);
             }
  frac_simplifie(X);
}

References B, D, frac::den, frac_simplifie(), FWD_OFL_CTRL, frac::num, tag, value_mult, value_protected_mult, value_substract, and X.

Referenced by full_pivot().

Here is the call graph for this function:

Here is the caller graph for this function:

full_pivot_sioux()

void full_pivot_sioux	(	frac *	X,
		frac	A,
		frac	B,
		frac	C,
		frac	D,
		bool	ofl_ctrl
	)

u*v*w == B*C/D

u*=v

u*=w

u*=v

u*=w

Definition at line 652 of file sc_simplex_feasibility_fixprec.c.

{
  frac u,v,w;
  tag("FULL_PIVOT_SIOUX")
    AFF(u,B); AFF(v,C); INV_ZERO_CTRL(w,D); /* u*v*w == B*C/D */
  frac_simpl(&(u.num),&(v.den));
  frac_simpl(&(u.num),&(w.den));
  frac_simpl(&(v.num),&(u.den));
  frac_simpl(&(v.num),&(w.den));
  frac_simpl(&(w.num),&(u.den));
  frac_simpl(&(w.num),&(v.den));
  if (ofl_ctrl == FWD_OFL_CTRL) {
    u.num = value_protected_mult(u.num,v.num); /* u*=v */
    u.den = value_protected_mult(u.den,v.den);
    u.num =value_protected_mult(u.num,w.num); /* u*=w */
    u.den =value_protected_mult(u.den,w.den);
  }
  else {
    u.num = value_mult(u.num,v.num); /* u*=v */
    u.den = value_mult(u.den,v.den);
    u.num =value_mult(u.num,w.num); /* u*=w */
    u.den =value_mult(u.den,w.den);
  }
  frac_sub(X,A,u,ofl_ctrl);
 
}

References AFF, B, D, frac::den, frac_simpl(), frac_sub(), FWD_OFL_CTRL, INV_ZERO_CTRL, frac::num, tag, value_mult, value_protected_mult, and X.

Referenced by full_pivot().

Here is the call graph for this function:

Here is the caller graph for this function:

hash()

static int hash ( Variable  s )

static

dump_tableau

calcule le hashcode d'un pointeur sous forme d'un nombre compris entre 0 et MAX_VAR

Definition at line 867 of file sc_simplex_feasibility_fixprec.c.

{ long l ;
  l=(long)s % MAX_VAR ;
  return l ;
}

References MAX_VAR.

Referenced by sc_min(), and sc_simplex_feasibility_ofl_ctrl_fixprec().

Here is the caller graph for this function:

partial_pivot()

void partial_pivot	(	frac *	X,
		frac	B,
		frac	C,
		frac	D,
		bool	ofl_ctrl
	)

Definition at line 755 of file sc_simplex_feasibility_fixprec.c.

{
 
  if (direct_p(B.den) && direct_p(C.den) && direct_p(D.num))
    {
      partial_pivot_direct(X,B,C,D,ofl_ctrl);
    }
  else
    {
      partial_pivot_sioux(X,B,C,D,ofl_ctrl);
    }
}

References B, D, frac::den, direct_p, partial_pivot_direct(), partial_pivot_sioux(), and X.

Referenced by pivot().

Here is the call graph for this function:

Here is the caller graph for this function:

partial_pivot_direct()

void partial_pivot_direct	(	frac *	X,
		frac	B,
		frac	C,
		frac	D,
		bool	ofl_ctrl
	)

Definition at line 736 of file sc_simplex_feasibility_fixprec.c.

{ tag("PARTIAL_PIVOT_DIRECT") 
    if (ofl_ctrl == FWD_OFL_CTRL) {
      X->num = value_protected_mult(B.num,C.num);
      X->num = value_protected_mult(X->num,D.den);
      value_oppose(X->num);
      X->den =value_protected_mult(B.den,C.den);
      X->den = value_protected_mult(X->den,D.num);
    }
    else {
      X->num = value_mult(B.num,C.num);
      X->num = value_mult(X->num,D.den);
      value_oppose(X->num);
      X->den =value_mult(B.den,C.den);
      X->den = value_mult(X->den,D.num);
    }
  frac_simplifie(X);
}

References B, D, frac::den, frac_simplifie(), FWD_OFL_CTRL, frac::num, tag, value_mult, value_oppose, value_protected_mult, and X.

Referenced by partial_pivot().

Here is the call graph for this function:

Here is the caller graph for this function:

partial_pivot_sioux()

void partial_pivot_sioux	(	frac *	X,
		frac	B,
		frac	C,
		frac	D,
		bool	ofl_ctrl
	)

idem if A==0

u=simplify(b*c)

x=simplify(u/d)

x=-x

Definition at line 726 of file sc_simplex_feasibility_fixprec.c.

{
  tag("PARTIAL_PIVOT_SIOUX")
    frac u =  {(Value)0,(Value)1,0};
 
  frac_mul(&u,B,C,ofl_ctrl); /* u=simplify(b*c) */
  frac_div(X,u,D,ofl_ctrl); /* x=simplify(u/d) */
  value_oppose(X->num);   /* x=-x */
}

References B, D, frac_div(), frac_mul(), tag, value_oppose, and X.

Referenced by partial_pivot().

Here is the call graph for this function:

Here is the caller graph for this function:

pivot()

void pivot	(	frac *	X,
		frac	A,
		frac	B,
		frac	C,
		frac	D,
		bool	ofl_ctrl
	)

a==0?

b*c/d != 0, calculons!

a!=0

b*c/d != 0, calculons!

no den to compute

well, we must compute the full formula!

Definition at line 770 of file sc_simplex_feasibility_fixprec.c.

{
  if (value_zero_p(D.num))
    fprintf(stderr,"division of zero!!!");
  DEBUG3(fprintf(stdout, "pivot on: ");
         printfrac(A);
         printfrac(B);
         printfrac(C);
         printfrac(D));
  if (value_zero_p(A.num))/* a==0? */
    {
      if (value_zero_p(B.num) || value_zero_p(C.num) || value_zero_p(D.den))
        { MET_ZERO((*X)); }
      else /* b*c/d != 0, calculons! */
        { partial_pivot(X,B,C,D,ofl_ctrl);}
    }
  else /* a!=0 */
    if (value_zero_p(B.num) || value_zero_p(C.num) || value_zero_p(D.den))
      { AFF_PX(X,A);}
    else /*  b*c/d != 0, calculons! */
      if (value_one_p(D.num) && value_one_p(A.den) &&
          value_one_p(B.den) && value_one_p(C.den))
        { /* no den to compute */
          Value v;
          if (ofl_ctrl == FWD_OFL_CTRL) {
            v= value_protected_mult(B.num,C.num);
            v = value_protected_mult(v,D.den);
          }
          else {
            v= value_mult(B.num,C.num);
            v = value_mult(v,D.den);
          }
          X->num=value_minus(A.num,v);
          X->den=VALUE_ONE;
        }
      else /* well, we must compute the full formula! */
        { full_pivot(X,A,B,C,D,ofl_ctrl);}
  DEBUG3(fprintf(stdout, " = ");
         printfrac(X); fprintf(stdout, "\n"));
}

References AFF_PX, B, D, DEBUG3, frac::den, fprintf(), full_pivot(), FWD_OFL_CTRL, MET_ZERO, frac::num, partial_pivot(), printfrac(), value_minus, value_mult, VALUE_ONE, value_one_p, value_protected_mult, value_zero_p, and X.

Referenced by choisir_piv(), dualentier(), gauss(), iprimal(), is2(), pivoter(), sc_min(), sc_simplex_feasibility_ofl_ctrl_fixprec(), simbad(), simplexe(), and simplexedual().

Here is the call graph for this function:

Here is the caller graph for this function:

printfrac()

static void printfrac ( frac  x )

static

Definition at line 832 of file sc_simplex_feasibility_fixprec.c.

                              {
    printf(" "); print_Value(x.num);
    printf("/"); print_Value(x.den);
}

References print_Value(), printf(), and x.

Referenced by pivot(), sc_min(), and sc_simplex_feasibility_ofl_ctrl_fixprec().

Here is the call graph for this function:

Here is the caller graph for this function:

sc_simplex_feasibility_ofl_ctrl_fixprec()

bool sc_simplex_feasibility_ofl_ctrl_fixprec	(	Psysteme	sc,
		int	ofl_ctrl
	)

fonction de calcul de la faisabilite' d'un systeme d'equations et d'inequations Auteur : Robert Mahl, Date : janvier 1994

Retourne : 1 si le systeme est soluble (faisable) en rationnels, 0 s'il n'y a pas de solution.

overflow control : ofl_ctrl == NO_OFL_CTRL => no overflow control ofl_ctrl == FWD_OFL_CTRL
=> overflow control is made THROW(overflow_error,5) BC, 13/12/94

All the folowing automatic variables are used when coming back from longjmp (i.e. in a CATCH block) so they need to be declared volatile as specified by the documentation

tete de liste des noms de variables

Necessaire de declarer "hashtable" static pour initialiser tout automatiquement a` 0. Necessaire de chainer les enregistrements pour reinitialiser a 0 en sortie de la procedure.

tableau des inegalite's

les colonnes 0 et 1 sont reservees au terme const:

valeur retournee par feasible

objectif de max pour simplex : somme des (b2,c2) termes constants "inferieurs"

count number of calls of this function (at the beginning)

int soluble = 1;

N.

recompute the base so as to only allocate necessary columns some bases are quite large although all variables do not appear in actual contraints. The base is used to store all variants in preconditions for instance.

Allocation a priori du tableau des egalites. "eg" : tableau a "nb_eq" lignes et "dimension"+2 colonnes.

DEBUG(fprintf(stdout, "\n\n IN sc_simplexe_feasibility_ofl_ctrl:\n"); sc_fprint(stdout, sc, default_variable_to_string);)

c_default_dump_to_file(); print to file

the input Psysteme must be consistent; this is not the best way to do this; array bound checks should be added instead in proper places; no time to do it properly for the moment. BC.

Do not allocate place for NULL constraints

ifscdebug(2) { fprintf(stderr,"[sc_simplexe_feasibility_ofl_ctrl] arithmetic error\n"); }

restore initial base

ethrow whatever the exception is

HROW(user_exception_error);

need CATCH(user_exception_error) before calling sc_simplexe_feasibility)

if don't catch exception, then default is feasible

f (ofl_ctrl == FWD_OFL_CTRL)

HROW(overflow_error);

eturn true; default is feasible

of CATCH(simplex_arithmetic_error)

egin of TRY

Allocation a priori du tableau du simplex "t" par colonnes. Soit "dimension" le nombre de variables, la taille maximum du tableau est de (1 + nb_ineq) lignes et de (2 + dimension + nb_ineq + nb_eq) colonnes On y ajoute en fait le double du nombre d'egalite's. Ce tableau sera rempli de la facon suivante :

• ligne 0 : critere d'optimisation
• lignes 1 a nb_ineq : les inequations
• colonne 0 : le terme constant (composante de poids 1)
• colonne 1 : le terme constant (composante de poids M)
• colonnes 2 et suivantes : les elements initiaux et les termes d'ecart Le tableau a une derniere colonne temporaire pour pivoter un vecteur unitaire.

Initialisation de l'objectif

Entree des inegalites dans la table

skip if empty

val terme const

Creation de variable d'ecart ?

Mise a jour des colonnes 0 et 1

Element de poids M en 1ere colonne

Creation d'une colonne cachee

NON IMPLEMENTE'

Elimination de Gauss-Jordan dans le tableau "eg" Chaque variable a` eliminer est marquee eg[ ].existe = 2 Si le processus d'elimination ne revele pas d'impossibilite', il est suivi du processus d'elimination dans les inegalites.

FIN DE NON IMPLEMENTE'

SOLUTION PROVISOIRE Pour chaque egalite on introduit 2 inequations complementaires

Added by bc: do nothing for nul equalities

val terme const

Creation de variable d'ecart ?

Mise a jour des colonnes 0 et 1

Element de poids M en 1ere colonne

Creation d'une colonne cachee

Added by bc: do nothing for nul equalities

val terme const

Creation de variable d'ecart ?

Mise a jour des colonnes 0 et 1

Element de poids M en 1ere colonne

Creation d'une colonne cachee

FIN DE SOLUTION PROVISOIRE

Algorithme du simplexe - methode primale simple. L'objectif est d'etudier la faisabilite' d'un systeme de contraintes sans trouver l'optimum. Les contraintes ont la forme : a x <= b et d x = e La methode de resolution procede comme suit :

1) Creer un tableau a b d e Eliminer autant de variables que posible par Gauss-Jordan

2) Travailler sur les inegalites seulement. Poser x = x' - M 1 ou` 1 est le vecteur de chiffres 1. Les inequations prennent alors la forme : a1 x <= b1 + M c1 a2 x >= b2 + M c2 avec c1 et c2 positifs On introduit les variables d'ecart y (autant que d'inequations du 2eme type) et on cherche max{1(a2 x - y) | x,y >= 0 ; a1 x <= b1 + M c1 ; a2 x - y <= b2 + M c2} On cree donc le tableau : 0 0 1 a2 1 0 0 b1 c1 a1 0 I 0 b2 c2 a2 -I 0 I

On applique ensuite l'algorithme du simplex en se souvenant que c1 et c2 sont a multiplier par un infiniment grand. Si l'optimum est egal a (1 b2 , 1 c2), il y a une solution.

Structures de donnees : on travaille sur des tableaux de fractions rationnelles.

Recherche d'un nombre negatif 1ere ligne

Terminaison

Recherche de la ligne de pivot
si ii==-1, pas encore trouve, min{1,2} non valides...

computing rapport{1,2}

first assignment is forced

POSITIF(t[jj].colonne[i]))

it may skip over

Cas d'impossibilite'

Modification des numeros des variables

Pivot autour de la ligne ii / colonne jj Dans ce qui suit, j = colonne courante, k = numero element dans la nouvelle colonne qui remplacera la colonne j, cc = element (colonne j, ligne ii)

Remplir la derniere colonne temporaire de t qui contient un 1 en position ligne ii

N

0 en colonne j ligne t[jj].colonne[i1].numero

0 en col jj, ligne t[j].colonne[i].numero

Restauration des entrees vides de la table hashee

restore initial base

Definition at line 888 of file sc_simplex_feasibility_fixprec.c.

{
    Pcontrainte pc, pc_tmp ;
    Pvecteur pv ;
    /* All the folowing automatic variables are used when coming back from
     * longjmp (i.e. in a CATCH block) so they need to be declared volatile as
     * specified by the documentation*/      
    intptr_t volatile premier_hash = PTR_NIL; /* tete de liste des noms de variables */
    /* Necessaire de declarer "hashtable" static 
     *  pour initialiser tout automatiquement a` 0.
     * Necessaire de chainer les enregistrements
     *  pour reinitialiser a 0
     *  en sortie de la procedure.
     */
    static hashtable_t volatile hashtable[MAX_VAR];
    Pbase volatile saved_base;
    int volatile saved_dimension;
    // tableau * volatile eg = NULL; /* tableau des egalite's  */
    tableau * volatile t = NULL; /* tableau des inegalite's  */
    frac * volatile nlle_colonne = NULL;
    /* les colonnes 0 et 1 sont reservees au terme const: */
    int compteur = 2 ;
    intptr_t i, j, k, h, trouve, ligne, i0, i1, jj, ii ;
    Value poidsM, valeur, tmpval;
    intptr_t w ;
    int soluble; /* valeur retournee par feasible */
    frac* colo;
    frac objectif[2] ; /* objectif de max pour simplex : 
                          somme des (b2,c2) termes constants "inferieurs" */
    frac rapport1, rapport2, min1, min2, piv, cc ;
 
    
    DEBUG(static int simplex_sc_counter = 0;)
      /* count number of calls of this function (at the beginning)       */
 
      soluble=1;/* int soluble = 1;*/
 
    rapport1 =frac0, rapport2 =frac0, min1 =frac0, min2 =frac0, piv=frac0, cc =frac0 ;
    objectif[0] = frac0, objectif[1] = frac0;
    i=-1, j=-1, k=-1, h=-1, trouve=-1, ligne=-1, i0=-1, i1=-1, jj=-1, ii=-1;
    poidsM =-1, valeur=-1, tmpval=-1,w=-1;/*DN.*/
 
    /* recompute the base so as to only allocate necessary columns
     * some bases are quite large although all variables do not appear in
     * actual contraints. The base is used to store all variants in
     * preconditions for instance.
     */
  
    saved_base = sc_base(sc);
    saved_dimension = sc_dimension(sc);
    sc_base(sc) = BASE_NULLE;
 
    sc_creer_base(sc);
 
 
    /* Allocation a priori du tableau des egalites.
     * "eg" : tableau a "nb_eq" lignes et "dimension"+2 colonnes.
     */
    if (ofl_ctrl!=FWD_OFL_CTRL)
        fprintf(stderr, "[sc_simplexe_feasibility_ofl_ctrl] "
                "should not (yet) be called with control %d...\n", ofl_ctrl);
 
    /* DEBUG(fprintf(stdout, "\n\n IN sc_simplexe_feasibility_ofl_ctrl:\n");
       sc_fprint(stdout, sc, default_variable_to_string);)*/
 
    DEBUG(simplex_sc_counter ++;
                  fprintf(stderr,"BEGIN SIMPLEX : %d th\n",simplex_sc_counter);
                  sc_default_dump(sc);/*sc_default_dump_to_file(); print to file          */
          );
 
    /* the input Psysteme must be consistent; this is not the best way to
     * do this; array bound checks should be added instead in proper places;
     * no time to do it properly for the moment. BC.
     */
        linear_assert("sc is weakly consistent", sc_weak_consistent_p(sc));
 
    /* Do not allocate place for NULL constraints */
    NB_EQ = 0;
    NB_INEQ = 0;
    for(pc_tmp = sc->egalites; pc_tmp!= NULL; pc_tmp=pc_tmp->succ)
    {
        if (pc_tmp->vecteur != NULL)
            NB_EQ++;
    }
    for(pc_tmp = sc->inegalites; pc_tmp!= NULL; pc_tmp=pc_tmp->succ)
    {
        if (pc_tmp->vecteur != NULL)
            NB_INEQ++;
    }
 
    CATCH(simplex_arithmetic_error|timeout_error|overflow_error)
    {
      /*      ifscdebug(2) {
        fprintf(stderr,"[sc_simplexe_feasibility_ofl_ctrl] arithmetic error\n");
        }
      */
      DEBUG(fprintf(stderr, "arithmetic error or timeout in simplex\n"););
 
        for(i = premier_hash ; i != PTR_NIL; i = hashtable[i].succ) {
          hashtable[i].nom =  VARIABLE_UNDEFINED ;
          hashtable[i].numero = 0 ;
          hashtable[i].hash = 0 ;
          hashtable[i].val = (Value)0 ;
        }
 
      for(i=0;i<(3 + NB_INEQ + NB_EQ + DIMENSION); i++)
        free(t[i].colonne);
      free(t);
      free(nlle_colonne);
 
      /* restore initial base */
      base_rm(sc_base(sc));
      sc_base(sc) = saved_base;
      sc_dimension(sc) = saved_dimension;
 
#ifdef CONTROLING
      alarm(0); /*clear the alarm*/
#endif
 
      if (ofl_ctrl == FWD_OFL_CTRL) {
        RETHROW(); /*rethrow whatever the exception is*/
      }
      /*THROW(user_exception_error);*/
      /* need CATCH(user_exception_error) before calling sc_simplexe_feasibility)*/
      ifscdebug(5) {fprintf(stderr,"DNDNDN WARNING: Exception not treated, return feasible!");}
      return true; /* if don't catch exception, then default is feasible */
 
      /*if (ofl_ctrl == FWD_OFL_CTRL)  */
      /*THROW(overflow_error);*/
 
      /*return true;  default is feasible */
    }/* of CATCH(simplex_arithmetic_error)*/
 
    /*begin of TRY*/
 
#ifdef CONTROLING
    /*start the alarm*/
    if (CONTROLING_TIMEOUT_SIMPLEX) {
      signal(SIGALRM, controling_catch_alarm_Simplex);
      alarm(CONTROLING_TIMEOUT_SIMPLEX);
    } /*else nothing*/
#endif
 
 
 
    /* Allocation a priori du tableau du simplex "t" par
     * colonnes. Soit
     * "dimension" le nombre de variables, la taille maximum
     * du tableau est de (1 + nb_ineq) lignes
     * et de (2 + dimension + nb_ineq + nb_eq) colonnes
     * On y ajoute en fait le double du nombre d'egalite's.
     * Ce tableau sera rempli de la facon suivante :
     * - ligne 0 : critere d'optimisation
     * - lignes 1 a nb_ineq : les inequations
     * - colonne 0 : le terme constant (composante de poids 1)
     * - colonne 1 : le terme constant (composante de poids M)
     * - colonnes 2 et suivantes : les elements initiaux
     *   et les termes d'ecart
     * Le tableau a une derniere colonne temporaire pour
     *  pivoter un vecteur unitaire.
     *     */
 
    t = (tableau*)malloc((3 + NB_INEQ + NB_EQ + DIMENSION)*sizeof(tableau));
    for(i=0;i<(3 + NB_INEQ + NB_EQ + DIMENSION); i++) {
        t[i].colonne= (frac*) malloc((4 + 2*NB_EQ + NB_INEQ)*sizeof(frac)) ;
        frac_init(t[i].colonne,4 + 2*NB_EQ + NB_INEQ);
        t[i].existe = 0 ;
        t[i].taille = 1 ;
        t[i].colonne[0].numero = 0 ;
        t[i].colonne[0].num = VALUE_ZERO ;
        t[i].colonne[0].den = VALUE_ONE ;
    }
    nbvariables= 0 ;
    /* Initialisation de l'objectif */
 
    for(i=0;i<=1;i++)
        objectif[i].num=VALUE_ZERO, objectif[i].den=VALUE_ONE;
 
    DEBUG2(dump_hashtable(hashtable);)
 
    /* Entree des inegalites dans la table */
 
    for(pc=sc->inegalites, ligne=1; pc!=0; pc=pc->succ, ligne++)
    {
        pv=pc->vecteur;
        if (pv!=NULL) /* skip if empty */
        {
            valeur = VALUE_ZERO ;
            poidsM = VALUE_ZERO ;
            for(; pv !=0 ; pv=pv->succ)
                if(vect_coeff(pv->var,sc_base(sc)))
                    value_addto(poidsM,pv->val) ;
                else
                    valeur = value_uminus(pv->val) ; /* val terme const */
 
            for(pv=pc->vecteur ; pv !=0 ; pv=pv->succ) {
                if(value_notzero_p(vect_coeff(pv->var,sc_base(sc)))) {
                    h = hash((Variable)  pv->var) ; trouve=0 ;
                    while (VARIABLE_DEFINED_P(hashtable[h].nom))  {
                        if (hashtable[h].nom==pv->var) {
                            trouve=1 ;
                            break ;
                        }
                        else { h = (h+1) % MAX_VAR ; }
                    }
                    if(!trouve) {
                        hashtable[h].succ=premier_hash ;
                        premier_hash = h ;
                        hashtable[h].val = VALUE_NAN ;
                        hashtable[h].numero=compteur++ ;
                        hashtable[h].nom=pv->var ;
                        CREVARVISIBLE ;
                    }
                    linear_assert("current NUMERO in bound",
                                                  (NUMERO) < (3 + NB_INEQ + NB_EQ + DIMENSION));
                    if (value_neg_p(poidsM) || 
                        (value_zero_p(poidsM) && value_neg_p(valeur)))
                        {value_addto(t[NUMERO].colonne[0].num,pv->val),
                           t[NUMERO].colonne[0].den = VALUE_ONE ;}
                    t[NUMERO].existe = 1 ;
                    t[NUMERO].colonne[t[NUMERO].taille].numero=ligne ;
                    if(value_neg_p(poidsM) || 
                       (value_zero_p(poidsM) && value_neg_p(valeur)))
                        tmpval = value_uminus(pv->val) ; 
                    else tmpval = pv->val ;
                    t[NUMERO].colonne[t[NUMERO].taille].num = tmpval ;
                    t[NUMERO].colonne[t[NUMERO].taille].den = VALUE_ONE ;
                    t[NUMERO].taille++ ;
                }
            }
 
            /* Creation de variable d'ecart ? */
            if(value_neg_p(poidsM) ||
               (value_zero_p(poidsM) && value_neg_p(valeur))) {
                DEBUG1(dump_tableau("cre var ec", t, compteur);)
                i=compteur++ ;
                CREVARVISIBLE ;
                t[i].existe = 1 ; t[i].taille = 2 ;
                t[i].colonne[0].num = VALUE_ONE ;
                t[i].colonne[0].den = VALUE_ONE ;
                DEBUG1(printf("ligne ecart = %ld, colonne %ld\n",ligne,i);)
                t[i].colonne[1].numero = ligne ;
                t[i].colonne[1].num = VALUE_MONE ;
                t[i].colonne[1].den = VALUE_ONE ;
                value_oppose(poidsM);
                value_oppose(valeur);
                value_addto(objectif[0].num,valeur) ; 
                value_addto(objectif[1].num,poidsM) ;
            }
 
            /* Mise a jour des colonnes 0 et 1 */
            t[0].colonne[t[0].taille].numero = ligne ;
            t[0].colonne[t[0].taille].den = VALUE_ONE ;
            t[0].colonne[t[0].taille].num = valeur ;
            t[0].existe = 1 ;
            t[0].taille++ ;
            /* Element de poids M en 1ere colonne */
            t[1].colonne[t[1].taille].numero = ligne ;
            t[1].colonne[t[1].taille].num = poidsM ;
            t[1].colonne[t[1].taille].den = VALUE_ONE ;
            t[1].existe = 1 ;
            t[1].taille++ ;
            /* Creation d'une colonne cachee */
            CREVARCACHEE ;
            DEBUG1(dump_tableau("cre col cach", t, compteur);)
                }
        else
            ligne--;
    }
 
    DEBUG1(dump_hashtable(hashtable);)
    DEBUG1(dump_tableau("avant sol prov", t, compteur);)
      
    /* NON IMPLEMENTE' */
    
    /* Elimination de Gauss-Jordan dans le tableau "eg"
     *  Chaque variable a` eliminer est marquee
     *  eg[ ].existe = 2
     *  Si le processus d'elimination ne revele pas
     *  d'impossibilite', il est suivi du processus
     *  d'elimination dans les inegalites.
     */
    /* FIN DE NON IMPLEMENTE' */
    
    /* SOLUTION PROVISOIRE
     *  Pour chaque egalite on introduit
     *  2 inequations complementaires
     */
    
    for(pc=sc->egalites ; pc!=0; pc=pc->succ, ligne++)
    {
        /* Added by bc: do nothing for nul equalities */
        if (pc->vecteur == NULL) continue;
 
        valeur = VALUE_ZERO ;
        poidsM = VALUE_ZERO ;
        for(pv=pc->vecteur ; pv !=0 ; pv=pv->succ)
            if(vect_coeff(pv->var,sc_base(sc)))
                value_addto(poidsM,pv->val) ;
            else valeur = value_uminus(pv->val); /* val terme const */
        for(pv=pc->vecteur ; pv !=0 ; pv=pv->succ) {
            if(vect_coeff(pv->var,sc_base(sc))) {
                h = hash((Variable) pv->var) ; trouve=0 ;
                while (VARIABLE_DEFINED_P(hashtable[h].nom))  {
                    if (hashtable[h].nom==pv->var) {
                        trouve=1 ;
                        break ;
                    }
                    else { h = (h+1) % MAX_VAR ; }
                }
                if(!trouve) {
                    hashtable[h].succ=premier_hash ;
                    premier_hash = h ;
                    hashtable[h].val = VALUE_NAN ;
                    hashtable[h].numero=compteur++ ;
                    CREVARVISIBLE ;
                    hashtable[h].nom=pv->var ;
                }
                                linear_assert("current NUMERO in bound",
                                  (NUMERO) < (3 + NB_INEQ + NB_EQ + DIMENSION));
                if(value_neg_p(poidsM) ||
                   (value_zero_p(poidsM) && value_neg_p(valeur)))
                    {value_addto(t[NUMERO].colonne[0].num,pv->val),
                       t[NUMERO].colonne[0].den = VALUE_ONE ;}
                t[NUMERO].existe = 1 ;
                t[NUMERO].colonne[t[NUMERO].taille].numero=ligne ;
                if(value_neg_p(poidsM) ||
                   (value_zero_p(poidsM) && value_neg_p(valeur)))
                    tmpval = value_uminus(pv->val);
                else tmpval = pv->val ;
                t[NUMERO].colonne[t[NUMERO].taille].num = tmpval ;
                t[NUMERO].colonne[t[NUMERO].taille].den = VALUE_ONE ;
                t[NUMERO].taille++ ;
            }
        }
        /* Creation de variable d'ecart ? */
        if(value_neg_p(poidsM) ||
           (value_zero_p(poidsM) && value_neg_p(valeur))) {
            i=compteur++ ;
            CREVARVISIBLE ;
            t[i].existe = 1 ; t[i].taille = 2 ;
            t[i].colonne[0].num = VALUE_ONE ;
            t[i].colonne[0].den = VALUE_ONE ;
            t[i].colonne[1].numero = ligne ;
            t[i].colonne[1].num = VALUE_MONE ;
            t[i].colonne[1].den = VALUE_ONE ;
            value_oppose(poidsM), 
            value_oppose(valeur);
            value_addto(objectif[0].num,valeur) ;
            value_addto(objectif[1].num,poidsM) ;
        }
        /* Mise a jour des colonnes 0 et 1 */
        t[0].colonne[t[0].taille].numero = ligne ;
        t[0].colonne[t[0].taille].num = valeur ;
        t[0].colonne[t[0].taille].den = VALUE_ONE ;
        t[0].existe = 1 ;
        t[0].taille++ ;
        /* Element de poids M en 1ere colonne */
        t[1].colonne[t[1].taille].numero = ligne ;
        t[1].colonne[t[1].taille].num = poidsM ;
        t[1].colonne[t[1].taille].den = VALUE_ONE ;
        t[1].existe = 1 ;
        t[1].taille++ ;
        /* Creation d'une colonne cachee */
        CREVARCACHEE ;
        DEBUG1(dump_tableau("cre col cach 2", t, compteur);)
    }
    
    for(pc=sc->egalites ; pc!=0; pc=pc->succ, ligne++)
    {
        /* Added by bc: do nothing for nul equalities */
        if (pc->vecteur == NULL) continue;
 
        valeur = VALUE_ZERO ;
        poidsM = VALUE_ZERO ;
        for(pv=pc->vecteur ; pv !=0 ; pv=pv->succ)
            if(vect_coeff(pv->var,sc_base(sc)))
                value_substract(poidsM, pv->val) ;
            else 
                valeur = pv->val ; /* val terme const */
 
        for(pv=pc->vecteur ; pv !=0 ; pv=pv->succ) {
            if (vect_coeff(pv->var,sc_base(sc))) {
                h = hash((Variable) pv->var) ; trouve=0 ;
                while (VARIABLE_DEFINED_P(hashtable[h].nom))  {
                    if (hashtable[h].nom==pv->var) {
                        trouve=1 ;
                        break ;
                    }
                    else { h = (h+1) % MAX_VAR ; }
                }
                if(!trouve) {
                    hashtable[h].succ=premier_hash ;
                    premier_hash = h ;
                    hashtable[h].val = VALUE_NAN ;
                    hashtable[h].numero=compteur++ ;
                    hashtable[h].nom=pv->var ;
                    CREVARVISIBLE ;
                }
                assert((NUMERO) < (3 + NB_INEQ + NB_EQ + DIMENSION));
                if(value_neg_p(poidsM) || 
                   (value_zero_p(poidsM) && value_neg_p(valeur)))
                    {value_substract(t[NUMERO].colonne[0].num,pv->val),
                       t[NUMERO].colonne[0].den = VALUE_ONE ;}
                t[NUMERO].existe = 1 ;
                t[NUMERO].colonne[t[NUMERO].taille].numero=ligne ;
                if(value_neg_p(poidsM) || 
                   (value_zero_p(poidsM) && value_neg_p(valeur)))
                    tmpval = pv->val ; 
                else tmpval = value_uminus(pv->val) ;
                t[NUMERO].colonne[t[NUMERO].taille].num = tmpval ;
                t[NUMERO].colonne[t[NUMERO].taille].den = VALUE_ONE ;
                t[NUMERO].taille++ ;
            }
        }
        /* Creation de variable d'ecart ? */
        if(value_neg_p(poidsM) || 
           (value_zero_p(poidsM) && value_neg_p(valeur))) {
            i=compteur++ ;
            CREVARVISIBLE ;
            t[i].existe = 1 ; t[i].taille = 2 ;
            t[i].colonne[0].num = VALUE_ONE ;
            t[i].colonne[0].den = VALUE_ONE ;
            t[i].colonne[1].numero = ligne ;
            t[i].colonne[1].num = VALUE_MONE ;
            t[i].colonne[1].den = VALUE_ONE ;
            value_oppose(poidsM), 
            value_oppose(valeur);
            value_addto(objectif[0].num,valeur) ;
            value_addto(objectif[1].num,poidsM) ;
        }
        /* Mise a jour des colonnes 0 et 1 */
        t[0].colonne[t[0].taille].numero = ligne ;
        t[0].colonne[t[0].taille].num = valeur ;
        t[0].colonne[t[0].taille].den = VALUE_ONE ;
        t[0].existe = 1 ;
        t[0].taille++ ;
        /* Element de poids M en 1ere colonne */
        t[1].colonne[t[1].taille].numero = ligne ;
        t[1].colonne[t[1].taille].num = poidsM ;
        t[1].colonne[t[1].taille].den = VALUE_ONE ;
        t[1].existe = 1 ;
        t[1].taille++ ;
        /* Creation d'une colonne cachee */
        CREVARCACHEE ;
        DEBUG1(dump_tableau("cre col cach 3", t, compteur);)
    }
    
    /* FIN DE SOLUTION PROVISOIRE */
    
    /* Algorithme du simplexe - methode primale simple.
     * L'objectif est d'etudier la faisabilite' d'un systeme
     * de contraintes sans trouver l'optimum.
     *   Les contraintes ont la forme : a x <= b
     *      et  d x = e
     * La methode de resolution procede comme suit :
     *
     *  1) Creer un tableau
     *       a  b
     *       d  e
     *     Eliminer autant de variables que posible par
     *    Gauss-Jordan
     *
     *  2) Travailler sur les inegalites seulement.
     *      Poser  x = x' - M 1
     *    ou` 1 est le vecteur de chiffres 1.
     *     Les inequations prennent alors la forme :
     *      a1 x <= b1 + M c1
     *      a2 x >= b2 + M c2
     *    avec c1 et c2 positifs
     *     On introduit les variables d'ecart y (autant que 
     *    d'inequations du 2eme type) et on cherche
     *      max{1(a2 x - y) | x,y >= 0 ; a1 x <= b1 + M c1 ;
     *                                a2 x - y <= b2 + M c2}
     *     On cree donc le tableau :
     *        0  0  1 a2     1  0  0
     *        b1 c1  a1      0  I  0
     *        b2 c2  a2     -I  0  I
     *
     *     On applique ensuite l'algorithme du simplex en
     *    se souvenant que c1 et c2 sont a multiplier par un
     *    infiniment grand.
     *     Si l'optimum est egal a (1 b2 , 1 c2), il y a une
     *    solution.
     *
     * Structures de donnees : on travaille sur des tableaux
     * de fractions rationnelles.
     */
    nlle_colonne=(frac *) malloc((4 + 2*NB_EQ + NB_INEQ+1)*sizeof(frac)) ;
    frac_init(nlle_colonne,4 + 2*NB_EQ + NB_INEQ);
    while(1) {
 
        /*  Recherche d'un nombre negatif 1ere ligne  
         */
        for(j=2, jj= -1 ;j<compteur;j++)
            if(t[j].existe && NEGATIF(t[j].colonne[0]))
            {
                jj=j ; break ;
            }
        
        /*  Terminaison  */
        if(jj == -1) { 
            bool cond;
            DEBUG1({
                printf ("solution :\n") ;
                dump_tableau("sol", t, compteur) ;
                printf("objectif : "); printfrac(objectif[0]) ; 
                        printfrac(objectif[1]) ; printf("\n") ;
            });
            
            if (ofl_ctrl == FWD_OFL_CTRL)
                cond = EGALOFL(objectif[0],t[0].colonne[0]) &&
                    EGALOFL(objectif[1],t[1].colonne[0]);
            else
                cond = EGAL(objectif[0],t[0].colonne[0]) &&
                    EGAL(objectif[1],t[1].colonne[0]);
            
            if(cond)
            {
                DEBUG1(printf("Systeme soluble (faisable) en rationnels\n");)
                SOLUBLE(1)
            }
            else 
            {
                DEBUG1(printf("Systeme insoluble (infaisable)\n");)
                SOLUBLE(0)
            }
            DEBUG1(printf("fin\n");)
        }
 
        DEBUG1(printf("1 : jj= %ld\n",jj);
              dump_tableau("avant ch pivot", t, compteur);)
 
        DEBUG1(min1.num = 32700; min1.den=1; min2=min1;)
        
        /*  Recherche de la ligne de pivot  
         *  si ii==-1, pas encore trouve, min{1,2} non valides...
         */
        for(i=1, i0=1, i1=1, ii=-1 ; i<t[jj].taille ; )
        {
            bool cond;
 
            DEBUG1(fprintf(stdout, "itering i{,0,1} = %ld %ld %ld\n", 
                          i, i0, i1);)
 
            if(((i0<t[0].taille && t[jj].colonne[i].numero <= 
                 t[0].colonne[i0].numero)  || i0>=t[0].taille)
               && ((i1<t[1].taille && t[jj].colonne[i].numero <=
                    t[1].colonne[i1].numero) || i1>=t[1].taille)) {
                if( POSITIF(t[jj].colonne[i])) {
                    /* computing rapport{1,2} 
                     */
                    frac f1 = 
                        (i0<t[0].taille &&
                         t[jj].colonne[i].numero==t[0].colonne[i0].numero)?
                             t[0].colonne[i0]:frac0;
                    frac f2 = t[jj].colonne[i];
                    frac f3 =
                        (i1<t[1].taille && 
                         t[jj].colonne[i].numero==t[1].colonne[i1].numero)?
                             t[1].colonne[i1]:frac0;
                        
                    frac_div(&rapport1,f1,f2,ofl_ctrl);
                    frac_div(&rapport2,f3,f2,ofl_ctrl);
        
    
                    DEBUG1(fprintf(stdout, "rapports:");
                          printfrac(rapport1);
                          printfrac(min1);
                          printfrac(rapport2);
                          printfrac(min2);
                          fprintf(stdout, "\nand cond: ");)
 
                    if (ii==-1) 
                        cond = true; /* first assignment is forced */
                    else
                        if (ofl_ctrl == FWD_OFL_CTRL)
                            cond = INFOFL(rapport2,min2) ||
                                (EGALOFL(rapport2,min2) && 
                                 INFOFL(rapport1,min1));
                        else
                            cond = INF(rapport2,min2) || 
                                (EGAL(rapport2,min2) && 
                                 INF(rapport1,min1));
                    
                    DEBUG1(fprintf(stdout, "%d\n", cond);)
 
                    if (cond) {
                        AFF(min1,rapport1) ;
                        AFF(min2,rapport2) ;
                        AFF(piv,t[jj].colonne[i]) ;
                        frac_simplifie(&piv);
                        ii=t[jj].colonne[i].numero ;
                    }
                } /* POSITIF(t[jj].colonne[i])) */
                i++ ;
            }
            else {
                if(i0<t[0].taille && /* it may skip over */
                   t[jj].colonne[i].numero> t[0].colonne[i0].numero) i0++ ;
                if(i1<t[1].taille && 
                   t[jj].colonne[i].numero > t[1].colonne[i1].numero) i1++ ;
            }
            
            DEBUG1(printf("i=%ld i0=%ld i1=%ld   %d %d %d\n",
                         i,i0,i1,
                         i<t[jj].taille? t[jj].colonne[i].numero: -1,
                         i0<t[0].taille? t[0].colonne[i0].numero: -1,
                         i1<t[1].taille? t[1].colonne[i1].numero: -1);)
        }
 
        /* Cas d'impossibilite'  */
        if(ii==-1) {
            DEBUG1(dump_tableau("sol infinie", t, compteur);
                   fprintf(stderr,"Solution infinie\n");)
            SOLUBLE(1)
        }
 
        /* Modification des numeros des variables */
 
        j = variables[jj-2];
        k = variablescachees[ii-1];
        variables[jj-2] = k;
        variablescachees[ii-1] = j;
 
        DEBUG2({
            printf("Visibles :");
            for(j=0;j<compteur-2;j++)
                printf(" %d",variables[j]);
            printf("\nCachees :");
            for(j=0;j<nbvariables;j++)
                printf(" %d",variablescachees[j]);
            printf("\n");
        });
 
        /*  Pivot autour de la ligne ii / colonne jj
         * Dans ce qui suit, j = colonne courante,
         *  k = numero element dans la nouvelle colonne
         *     qui remplacera la colonne j,
         *  cc = element (colonne j, ligne ii)
         */
        DEBUG(printf("Pivoter %ld %ld\n",ii,jj);)
        
        /* Remplir la derniere colonne temporaire de t
         *   qui contient un 1 en position ligne ii
         */
        t[compteur].taille = 2 ;
        t[compteur].colonne[0].numero = 0 ;
        t[compteur].colonne[0].num = VALUE_ZERO ;
        t[compteur].colonne[0].den = VALUE_ONE ;
        t[compteur].colonne[1].numero = ii;
        t[compteur].colonne[1].num = VALUE_ONE ;
        t[compteur].colonne[1].den = VALUE_ONE ;
        t[compteur].existe = 1 ;
 
        for(j=0 ; j<=compteur ; j=(j==(jj-1))?(j+2):(j+1)) {
            if(t[j].existe)
            {
                k=0 ;
                cc.num= VALUE_ZERO ; 
                cc.den= VALUE_ONE ;
                for(i=1;i<t[j].taille;i++)
                    if(t[j].colonne[i].numero==ii)
                    { AFF(cc,t[j].colonne[i]); break ; }
                    else if(t[j].colonne[i].numero>ii)
                    {cc.num= VALUE_ZERO ; cc.den=VALUE_ONE ; break ; }
                for(i=0,i1=0;i<t[j].taille || i1<t[jj].taille ;) {
 
                    DEBUG2(printf("k=%ld, j=%ld, i=%ld i1=%ld\n",k,j,i,i1);
                           printf("fractions: ");
                           printfrac(t[j].colonne[i]) ;
                           printfrac(t[jj].colonne[i1]) ;
                           if (value_zero_p(t[jj].colonne[i1].den))
                           printf("ATTENTION fraction 0/0 ");/*DN*/
                           printfrac(cc);
                           printfrac(piv);)
                    
                    if(i<t[j].taille &&  
                       i1<t[jj].taille && 
                       t[j].colonne[i].numero == t[jj].colonne[i1].numero) 
                    {   
                        if(t[j].colonne[i].numero == ii) {
                            AFF(nlle_colonne[k],t[j].colonne[i]);
                        } else {
                            frac *n = &nlle_colonne[k],
                                 *a = &t[j].colonne[i],
                                 *b = &t[jj].colonne[i1];
                             pivot(n, (*a), (*b), cc, piv,ofl_ctrl);
                        }
                        
                        if(i==0||nlle_colonne[k].num!=0) {
                            nlle_colonne[k].numero = t[j].colonne[i].numero ;
                            k++ ;
                        }
 
                        i++ ; i1++ ;
                    }
                    else
                        if(i>=t[j].taille || 
                           (i1<t[jj].taille && 
                            t[j].colonne[i].numero > t[jj].colonne[i1].numero))
                        {  
                            DEBUG1(
                                if (i<t[j].taille) 
                                {
                                    printf("t[j].colonne[i].numero > "
                                           "t[jj].colonne[i1].numero , "
                                           "k=%ld, j=%ld, i=%ld i1=%ld\n",
                                           k,j,i,i1);
                                    printf("j = %ld  t[j].taille=%d , "
                                           "t[jj].taille=%d\n",
                                           j,t[j].taille,t[jj].taille);
                                    printf("t[j].colonne[i].numero=%d , "
                                           "t[jj].colonne[i1].numero=%d\n",
                                           t[j].colonne[i].numero,
                                           t[jj].colonne[i1].numero);
                                });
 
                        /* 0 en colonne j  ligne t[jj].colonne[i1].numero */
                            if(t[jj].colonne[i1].numero == ii) {
                                AFF(nlle_colonne[k],frac0)
                            } else {
                                frac *n = &(nlle_colonne[k]),
                                     *b = &(t[jj].colonne[i1]);
                                pivot(n,frac0,(*b),cc,piv,ofl_ctrl) ;
                            }
 
                            if(i==0||nlle_colonne[k].num!=0)
                            {
                                nlle_colonne[k].numero = 
                                    t[jj].colonne[i1].numero ;
                                k++ ;
                            }
                            if(i1<t[jj].taille) i1++ ; else i++ ;
                        }
                        else if(i1>=t[jj].taille || 
                                t[j].colonne[i].numero < 
                                t[jj].colonne[i1].numero)
                        {
                            /* 0 en col jj, ligne t[j].colonne[i].numero */
                            DEBUG2(printf("t[j].colonne[i].numero < "
                                          "t[jj].colonne[i1].numero , "
                                          "k=%ld, j=%ld, i=%ld i1=%ld\n",
                                          k,j,i,i1);
                                   printf("j = %ld  t[j].taille=%d , "
                                          "t[jj].taille=%d\n",
                                          j,t[j].taille,t[jj].taille););
                            AFF(nlle_colonne[k],t[j].colonne[i]) ;
                            if(i==0||nlle_colonne[k].num!=0) {
                                nlle_colonne[k].numero = 
                                    t[j].colonne[i].numero ;
                                k++ ;
                            }
                            if(i<t[j].taille) i++ ; else i1++ ;
                        }
                        else
                            DEBUG2(printf(" ??? ");)
 
                    DEBUG2(printf(" -> ");
                           printfrac(nlle_colonne[k-1]);
                           printf(" [ligne numero %d]\n", 
                                  nlle_colonne[k-1].numero);)
 
                }
                if(j==compteur) w = jj ; else w = j ;
                colo = t[w].colonne ;
                t[w].colonne=nlle_colonne ;
                nlle_colonne = colo ;
                t[w].taille=k ;
                DEBUG1(printf("w = %ld  t[w].taille=%d \n",w,t[w].taille);
                      dump_tableau("last", t, compteur););
            }
        }
    }
 
    /* Restauration des entrees vides de la table hashee  */
    FINSIMPLEX :
    DEBUG1(dump_tableau("fin simplexe", t, compteur);)
    DEBUG(fprintf(stderr, "soluble = %d\n", soluble);)
 
    DEBUG(fprintf(stderr,"END SIMPLEX: %d th\n",simplex_sc_counter);)
 
    for(i = premier_hash ; i != PTR_NIL; i = hashtable[i].succ){
      hashtable[i].nom = VARIABLE_UNDEFINED;
      hashtable[i].numero = 0 ;
      hashtable[i].hash = 0 ;
      hashtable[i].val = (Value)0 ;
    }
 
    
    for(i=0;i<(3 + NB_INEQ + NB_EQ + DIMENSION); i++) 
      free(t[i].colonne);
    free(t);
    free(nlle_colonne);
 
#ifdef CONTROLING
    alarm(0); /*clear the alarm*/
#endif
    UNCATCH(simplex_arithmetic_error|timeout_error|overflow_error);
    
    /* restore initial base */
    vect_rm(sc_base(sc));
    sc_base(sc) = saved_base;
    sc_dimension(sc) = saved_dimension;
    
    return soluble;
}     /* main */

References AFF, assert, BASE_NULLE, base_rm, CATCH, col::colonne, CREVARCACHEE, CREVARVISIBLE, DEBUG, DEBUG1, DEBUG2, frac::den, DIMENSION, EGAL, Ssysteme::egalites, EGALOFL, col::existe, f2(), fprintf(), frac0, frac_div(), frac_init(), frac_simplifie(), free(), FWD_OFL_CTRL, hashtable_t::hash, hash(), Ssysteme::inegalites, INF, INFOFL, intptr_t, linear_assert, malloc(), MAX_VAR, NB_EQ, NB_INEQ, nbvariables, NEGATIF, hashtable_t::nom, frac::num, num, frac::numero, NUMERO, hashtable_t::numero, overflow_error, pivot(), POSITIF, printf(), printfrac(), PTR_NIL, RETHROW, sc_creer_base(), sc_default_dump(), sc_weak_consistent_p(), simplex_arithmetic_error, SOLUBLE, Scontrainte::succ, hashtable_t::succ, Svecteur::succ, col::taille, timeout_error, UNCATCH, hashtable_t::val, Svecteur::val, valeur(), value_addto, VALUE_MONE, VALUE_NAN, value_neg_p, value_notzero_p, VALUE_ONE, value_oppose, value_substract, value_uminus, VALUE_ZERO, value_zero_p, Svecteur::var, VARIABLE_DEFINED_P, VARIABLE_UNDEFINED, variables, variablescachees, vect_coeff(), vect_rm(), and Scontrainte::vecteur.

Referenced by sc_simplex_feasibility_ofl_ctrl().

Here is the call graph for this function:

Here is the caller graph for this function:

Variable Documentation

frac0

frac frac0 ={(Value)0,(Value)1,0}

static

Definition at line 508 of file sc_simplex_feasibility_fixprec.c.

Referenced by sc_min(), and sc_simplex_feasibility_ofl_ctrl_fixprec().

NB_EQ

int NB_EQ = 0

static

test du simplex : Si on compile grace a` "make simp" dans le repertoire /projects/C3/Linear/Development/polyedre/Tests alors on peut tester l'execution dans le meme directory en faisant : make test_simp

To replace #define NB_EQ and #define NB_INEQ - BC, 2/4/96 - NB_EQ and NB_INEQ are initialized at the beginning of the main subroutine. they represent the number of non-NULL constraints in sc. This is useful to allocate the minimum amount of memory necessary.

Definition at line 53 of file sc_simplex_feasibility_fixprec.c.

Referenced by sc_simplex_feasibility_ofl_ctrl_fixprec().

NB_INEQ

int NB_INEQ = 0

static

Definition at line 54 of file sc_simplex_feasibility_fixprec.c.

Referenced by sc_simplex_feasibility_ofl_ctrl_fixprec().

nbvariables

int nbvariables

static

Le nombre de variables visibles est : compteur-2 La i-eme variable visible a le numero : variables[i+1]=i (0 <= i < compteur-2) Le nombre de variables cachees est : nbvarables La i-eme variable cachee a le numero : variablescachees[i+1]=MAX_VAR+i-1 (0 <= i < nbvariables)

utilise'es par dump_tableau ; a rendre local

Definition at line 830 of file sc_simplex_feasibility_fixprec.c.

Referenced by sc_simplex_feasibility_ofl_ctrl_fixprec().

variables

int variables[MAX_VAR]

static

Definition at line 830 of file sc_simplex_feasibility_fixprec.c.

Referenced by build_and_test_dependence_context(), gfc2pips_namespace(), gfc2pips_vars_(), sc_restricted_to_variables_transitive_closure(), and sc_simplex_feasibility_ofl_ctrl_fixprec().

variablescachees

int variablescachees[MAX_VAR]

static

Definition at line 830 of file sc_simplex_feasibility_fixprec.c.

Referenced by sc_simplex_feasibility_ofl_ctrl_fixprec().