[Top][All Lists]
[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
[Help-glpk] Re: For complexity size isn't the issue (The case of the mis
From: |
Nigel Galloway |
Subject: |
[Help-glpk] Re: For complexity size isn't the issue (The case of the missing 12). |
Date: |
Fri, 4 Jul 2008 12:43:42 +0100 |
True, I could use some method to determine smallest multiple of
the largest common divisor, say Euclids Algorithm, but that is
rather what I want GLPK to do.
Note that if I want to do that then I just need to replace
a*x + b*y >= 1; with a*x + b*y = the gcd; The minimize is
then optional, and no cutting is required.
The range in which to minimize can be reduced:
param a, integer;
param b, integer;
var x, integer;
var y, integer;
s1: min(a,b)-1 >= a*x + b*y >= 1;
minimize s: a*x+b*y;
solve;
printf "%i %i\n",x,y;
data;
param a := 15;
param b := 35;
end;
--but the result is the same.
glpsol --math gcd.mathprog
Reading model section from gcd.mathprog...
Reading data section from gcd.mathprog...
18 lines were read
Generating s1...
Generating s...
Model has been successfully generated
glp_simplex: original LP has 2 rows, 2 columns, 4 non-zeros
Objective value = 1
OPTIMAL SOLUTION FOUND BY LP PRESOLVER
Integer optimization begins...
+ 0: mip = not found yet >= -inf (1; 0)
+ 5742: mip = not found yet >= 1.000000000e+00 (6456; 0)
+ 8489: mip = not found yet >= 1.000000000e+00 (9203; 0)
+ 10565: mip = not found yet >= 1.000000000e+00 (11279; 0)
+ 12305: mip = not found yet >= 1.000000000e+00 (13019; 0)
+ 13804: mip = not found yet >= 1.000000000e+00 (14518; 0)
+ 15133: mip = not found yet >= 1.000000000e+00 (15847; 0)
+ 16295: mip = not found yet >= 1.000000000e+00 (17009; 0)
+ 17319: mip = not found yet >= 1.000000000e+00 (18033; 0)
+ 18229: mip = not found yet >= 1.000000000e+00 (18943; 0)
+ 19058: mip = not found yet >= 1.000000000e+00 (19772; 0)
+ 19822: mip = not found yet >= 1.000000000e+00 (20536; 0)
+ 20538: mip = not found yet >= 1.000000000e+00 (21252; 0)
Time used: 60.0 secs. Memory used: 6.0 Mb.
+ 21211: mip = not found yet >= 1.000000000e+00 (21925; 0)
+ 21849: mip = not found yet >= 1.000000000e+00 (22563; 0)
.
.
.
+ 61982: mip = not found yet >= 9.999999999e-01 (62696; 0)
Time used: 840.0 secs. Memory used: 16.8 Mb.
+ 62148: mip = not found yet >= 9.999999999e-01 (62862; 0)
+ 62314: mip = not found yet >= 9.999999999e-01 (63028; 0)
+ 62479: mip = not found yet >= 9.999999999e-01 (63193; 0)
+ 62644: mip = not found yet >= 9.999999999e-01 (63358; 0)
+ 62809: mip = not found yet >= 9.999999999e-01 (63523; 0)
+ 62973: mip = not found yet >= 9.999999999e-01 (63687; 0)
+ 63137: mip = not found yet >= 9.999999999e-01 (63851; 0)
+ 63300: mip = not found yet >= 9.999999999e-01 (64014; 0)
+ 63463: mip = not found yet >= 9.999999999e-01 (64177; 0)
+ 63625: mip = not found yet >= 9.999999999e-01 (64339; 0)
+ 63786: mip = not found yet >= 9.999999999e-01 (64500; 0)
+ 63947: mip = not found yet >= 9.999999999e-01 (64661; 0)
Time used: 900.0 secs. Memory used: 17.3 Mb.
+ 64108: mip = not found yet >= 9.999999999e-01 (64822; 0)
+ 64268: mip = not found yet >= 9.999999999e-01 (64982; 0)
+ 64428: mip = not found yet >= 9.999999999e-01 (65142; 0)
+ 64587: mip = not found yet >= 9.999999999e-01 (65301; 0)
+ 64747: mip = not found yet >= 9.999999999e-01 (65461; 0)
+ 64905: mip = not found yet >= 9.999999999e-01 (65619; 0)
+ 65064: mip = not found yet >= 9.999999999e-01 (65778; 0)
+ 65221: mip = not found yet >= 9.999999999e-01 (65935; 0)
+ 65379: mip = not found yet >= 9.999999999e-01 (66093; 0)
+ 65536: mip = not found yet >= 9.999999999e-01 (66250; 0)
+ 65692: mip = not found yet >= 9.999999999e-01 (66406; 0)
+ 65848: mip = not found yet >= 9.999999999e-01 (66562; 0)
Time used: 960.0 secs. Memory used: 17.8 Mb.
+ 66004: mip = not found yet >= 9.999999999e-01 (66718; 0)
+ 66159: mip = not found yet >= 9.999999999e-01 (66873; 0)
^C
When the minimize is removed:
param a, integer;
param b, integer;
var x, integer;
var y, integer;
s1: min(a,b)-1 >= a*x + b*y >= 1;
/* minimize s: a*x+b*y; */
solve;
printf "%i %i\n",x,y;
data;
param a := 15;
param b := 35;
end;
Reading model section from gcd.mathprog...
Reading data section from gcd.mathprog...
18 lines were read
Generating s1...
Model has been successfully generated
glp_simplex: original LP has 1 row, 2 columns, 2 non-zeros
Objective value = 0
OPTIMAL SOLUTION FOUND BY LP PRESOLVER
Integer optimization begins...
+ 0: mip = not found yet >= -inf (1; 0)
+ 2: >>>>> 0.000000000e+00 >= 0.000000000e+00 0.0% (3; 0)
+ 2: mip = 0.000000000e+00 >= tree is empty 0.0% (0; 5)
INTEGER OPTIMAL SOLUTION FOUND
Time used: 0.0 secs
Memory used: 0.1 Mb (104365 bytes)
-2 1
Model has been successfully processed
which in this case is a correct solution, but isn't always. I used this
approach see in a C# program:
http://lists.gnu.org/archive/html/help-glpk/2008-02/msg00065.html
The same idea implemented in mathprog is:
param a, integer;
param b, integer;
param guess;
var x, integer;
var y, integer;
s1: 1 <= a*x+b*y <= guess;
solve;
param g := a*x + b*y;
param ga := a mod g;
param gb := b mod g;
param sol := if ga = 0 and gb = 0 then g else g-1;
printf "guess = %i\n",guess;
printf "x = %i, y = %i, a*x + b*y = %i\n",x,y,g;
printf "%s %i\n",if ga = 0 and gb = 0 then "Solution is " else "Not found try
",sol;
data;
param a := 2424;
param b := 772;
param guess := 771;
end;
glpsol --math axpby.lt.guess.mathprog
Reading model section from axpby.lt.guess.mathprog...
Reading data section from axpby.lt.guess.mathprog...
27 lines were read
Generating s1...
Model has been successfully generated
glp_simplex: original LP has 1 row, 2 columns, 2 non-zeros
Objective value = 0
OPTIMAL SOLUTION FOUND BY LP PRESOLVER
Integer optimization begins...
+ 0: mip = not found yet >= -inf (1; 0)
+ 4: >>>>> 0.000000000e+00 >= 0.000000000e+00 0.0% (4; 1)
+ 4: mip = 0.000000000e+00 >= tree is empty 0.0% (0; 9)
INTEGER OPTIMAL SOLUTION FOUND
Time used: 0.0 secs
Memory used: 0.1 Mb (104376 bytes)
guess = 771
x = 1, y = -3, a*x + b*y = 108
Not found try 107
Model has been successfully processed
Now if guess could be passed to glpsol as a parameter we could make the guess a
parameter (or if we could have recursive mathprog), but otherwise making the
changes to guess in a text editor we can proceed:
glpsol --math axpby.lt.guess.mathprog
Reading model section from axpby.lt.guess.mathprog...
Reading data section from axpby.lt.guess.mathprog...
27 lines were read
Generating s1...
Model has been successfully generated
glp_simplex: original LP has 1 row, 2 columns, 2 non-zeros
Objective value = 0
OPTIMAL SOLUTION FOUND BY LP PRESOLVER
Integer optimization begins...
+ 0: mip = not found yet >= -inf (1; 0)
+ 4: >>>>> 0.000000000e+00 >= 0.000000000e+00 0.0% (3; 3)
+ 4: mip = 0.000000000e+00 >= tree is empty 0.0% (0; 11)
INTEGER OPTIMAL SOLUTION FOUND
Time used: 0.0 secs
Memory used: 0.1 Mb (104376 bytes)
guess = 107
x = -7, y = 22, a*x + b*y = 16
Not found try 15
Model has been successfully processed
glpsol --math axpby.
lt.guess.mathprog
Reading model section from axpby.lt.guess.mathprog...
Reading data section from axpby.lt.guess.mathprog...
27 lines were read
Generating s1...
Model has been successfully generated
glp_simplex: original LP has 1 row, 2 columns, 2 non-zeros
Objective value = 0
OPTIMAL SOLUTION FOUND BY LP PRESOLVER
Integer optimization begins...
+ 0: mip = not found yet >= -inf (1; 0)
+ 14: >>>>> 0.000000000e+00 >= 0.000000000e+00 0.0% (5; 11)
+ 14: mip = 0.000000000e+00 >= tree is empty 0.0% (0; 31)
INTEGER OPTIMAL SOLUTION FOUND
Time used: 0.0 secs
Memory used: 0.1 Mb (104376 bytes)
guess = 15
x = -50, y = 157, a*x + b*y = 4
Solution is 4
Model has been successfully processed
This compares with Euclid's Algorithm where the sequence would be 772, 108, 16,
12, 4. So this program is not just Euclid's Algorith in disguise because GLPK
misses the 12.
Now if I try a guess of 3, there is no combination of x and y which produce a
solution but:
glpsol --nopresol --math axpby.lt.guess.mathprog
Reading model section from axpby.lt.guess.mathprog...
Reading data section from axpby.lt.guess.mathprog...
27 lines were read
Generating s1...
Model has been successfully generated
lpx_adv_basis: size of triangular part = 1
! 0: objval = 0.000000000e+00 infeas = 0.000000000e+00
OPTIMAL SOLUTION FOUND
Integer optimization begins...
+ 0: mip = not found yet >= -inf (1; 0)
+ 11168: mip = not found yet >= 0.000000000e+00 (5663; 0)
+ 16446: mip = not found yet >= 0.000000000e+00 (8302; 0)
+ 20422: mip = not found yet >= 0.000000000e+00 (10290; 0)
+ 23730: mip = not found yet >= 0.000000000e+00 (11944; 0)
+ 26624: mip = not found yet >= 0.000000000e+00 (13391; 0)
+ 29148: mip = not found yet >= 0.000000000e+00 (14653; 0)
+ 31314: mip = not found yet >= 0.000000000e+00 (15736; 0)
+ 33140: mip = not found yet >= 0.000000000e+00 (16649; 0)
+ 34734: mip = not found yet >= 0.000000000e+00 (17446; 0)
+ 36150: mip = not found yet >= 0.000000000e+00 (18154; 0)
^C
Which since there is no solution, optimal or not, I assume looking for it will
never stop.
> ----- Original Message -----
> From: address@hidden
> To: address@hidden
> Subject: Re: For complexity size isn't the issue.
> Date: Mon, 30 Jun 2008 15:11:22 -0700
>
>
> Hello Nigel!
>
> The search space is infinite with x, y in [-inf, +inf] with
> infinitely many global optimum solutions hence solution time for
> branch and bound will be infinite too.
>
> The optimum is the smallest multiple of
> the largest common divisor of a and b
> > = the right hand side of st1.
>
> Obviously the objective function will always be integer. Hence it
> would make sense to let the branch and bound add an extra cut
> limiting the objective function to the best found solution minus
> one.
>
> Together with a cut providing a correct lower limit this would let
> the search stop at the first optimum solution.
>
> Best regards
>
> Xypron
>
--
_______________________________________________
Surf the Web in a faster, safer and easier way:
Download Opera 9 at http://www.opera.com
Powered by Outblaze
[Prev in Thread] |
Current Thread |
[Next in Thread] |
- [Help-glpk] Re: For complexity size isn't the issue (The case of the missing 12).,
Nigel Galloway <=