Using simplex package to solve a linear program
- Subject: Using simplex package to solve a linear program
- From: Andrej Vodopivec
- Date: Sat, 14 May 2011 07:22:55 +0200
The implementation of the simplex algorithm is not very sophisticated.
You get an error, probably because of numerical instability.
You can convert the coefficients in the problem to rational numbers
and that will give you a better solution:
(%i1) load(simplex)$
(%i2) obj : 7.3*x1 + 6.9*x2 + 7.3*x3 + 7.5*x4 + 7.6*x5 + 6.0*x6 +
5.8*x7 + 4.3*x8 + 4.1*x9$
(%i3) eq1 : x1+x2+x3+x4+x5+x6+x7+x8+x9 = 1$
(%i4) eq2 : .2*x1 + .5*x2 + .3*x3 + .3*x4+ .3*x5 + .6*x6 + .4*x7 +
.1*x8 + .1*x9 = .3$
(%i5) eq3 : .3*x1 + .4*x2 + .2*x3 + .4*x4 + .3*x5 + .3*x6 + .5*x7 +
.3*x8 + .1*x9 = .3$
(%i6) eq4 : .5*x1 + .1*x2 + .5*x3 + .3*x4 + .4*x5 + .1*x6 + .1*x7 +
.6*x8 + .8*x9 = .4$
(%i7) eqs : [eq1,eq2,eq3,eq4]$
(%i8) [opt_val, sol]: minimize_lp( ratsimp(obj), ratsimp(eqs), all),
ratprint=false;
(%o8) [249/50,[x9=0,x8=3/5,x7=0,x6=2/5,x5=0,x4=0,x3=0,x2=0,x1=0]]
(%i9) %, numer;
(%o9) [4.98,[x9=0,x8=0.6,x7=0,x6=0.4,x5=0,x4=0,x3=0,x2=0,x1=0]]
(%i10) opt_val - obj, sol;
(%o10) 0.0
(%i11) eqs, sol;
(%o11) [1=1,0.3=0.3,0.3=0.3,0.4=0.4]
Lp, Andrej
On Fri, May 13, 2011 at 10:37 PM, Donald Bindner <dbindner at truman.edu> wrote:
> I was starting on Dantzig's undegraduate text in linear programming, and
> wanted to work through the first exercise in Maxima. ?I admit to being
> pretty new at this, but the Maxima documentation seemed pretty clear. ?Yet I
> get answers that I suspect to be incorrect. ?This is my session:
> (%i1) load( "simplex" )$
> (%i2) obj : 7.3*x1 + 6.9*x2 + 7.3*x3 + 7.5*x4 + 7.6*x5 + 6.0*x6 + 5.8*x7 +
> 4.3*x8 + 4.1*x9;
> (%i3) eq1 : x1+x2+x3+x4+x5+x6+x7+x8+x9 = 1;
> (%i4) eq2 : .2*x1 + .5*x2 + .3*x3 + .3*x4+ .3*x5 + .6*x6 + .4*x7 + .1*x8 +
> .1*x9 = .3;
> (%i5) eq3 : .3*x1 + .4*x2 + .2*x3 + .4*x4 + .3*x5 + .3*x6 + .5*x7 + .3*x8 +
> .1*x9 = .3;
> (%i6) eq4 : .5*x1 + .1*x2 + .5*x3 + .3*x4 + .4*x5 + .1*x6 + .1*x7 + .6*x8 +
> .8*x9 = .4;
> (%i7) eqs : [eq1,eq2,eq3,eq4];
> (%i8) minimize_lp( obj, eqs), nonegative_lp=true;
> (%o8) [4.95, [x9 = - 0.10714285714286, x8 = 0.75, x7 = - 0.10714285714286,
> ? ? ? ? ? ? ? ? x6 = 0.46428571428571, x5 = 0, x4 = 0, x3 = 0, x2 = 0, x1 =
> 0]]
> According to the documentation, nonegative_lp=true is supposed to assume
> that all decision variables are positive, yet in this solution x9 and x7 are
> clearly negative.
> If we instead specify the constraints by hand we get:
> (%i9)?eqs2 : [eq1,eq2,eq3,eq4,x1>0,x2>0,x3>0,x4>0,x5>0,x6>0,x7>0,x8>0,x9>0];
> (%i10) minimize_lp(obj, eqs2);
> (%o10) [7.6, [x9 = 0.0, x8 = 0.0, x7 = 0, x6 = 0, x5 = 1.0, x4 = 0.0,
> ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?x3 = 0.0, x2 = 0.0, x1 =
> 0.0]]
> This solution has all non-negative variables, however it isn't really a good
> minimum. ?We can verify, for example that there are better solutions:
> (%i11) obj, x6=.25,x7=.25,x8=.25,x9=.25, x1=0,x2=0,x3=0,x4=0,x5=0;
> (%o11) ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? 5.05
> (%i12) eqs, x6=.25,x7=.25,x8=.25,x9=.25, x1=0,x2=0,x3=0,x4=0,x5=0;
> (%o12) ? ? ? ? ? ?[1.0 = 1, 0.3 = 0.3, 0.3 = 0.3, 0.4 = 0.4]
> Am I doing something wrong, or is the simplex package not working?
> I have verified this behavior in both Maxima 5.22post on an Ubuntu Maverick
> system, and in Maxima?Maxima 5.20.1 on an Ubuntu Lucid system.
> Thanks,
> Don
>
>
>
> _______________________________________________
> Maxima mailing list
> Maxima at math.utexas.edu
> http://www.math.utexas.edu/mailman/listinfo/maxima
>
>