Opened 6 years ago
Last modified 2 months ago
#18735 needs_review enhancement
MixedIntegerLinearProgram/HybridBackend: Reconstruct exact rational/algebraic basic solution
Reported by:  mkoeppe  Owned by:  

Priority:  major  Milestone:  sage9.5 
Component:  numerical  Keywords:  lp 
Cc:  yzh, ncohen, dimpase  Merged in:  
Authors:  Matthias Koeppe, Yuan Zhou  Reviewers:  
Report Upstream:  N/A  Work issues:  
Branch:  u/yzh/hybrid_backend (Commits, GitHub, GitLab)  Commit:  71ae94d57b06bd58a1cf02a9b286bcbe881a0fcc 
Dependencies:  #18685, #18688, #20296  Stopgaps: 
Description (last modified by )
Sometimes one can use a fast numerical LP solver to solve a problem to "optimality",
then reconstruct the primal and dual solution in rational arithmetic (or over whatever base_ring was used...) and in this way prove that this basis is indeed optimal.
MixedIntegerLinearProgram
should support this mode of operation.
The current branch, on top of #20296, attempts to do this by implementing a HybridBackend
, which delegates to two backends:
 a fast, possibly inexact backend (Gurobi or GLPK or even GLPK with glp_exact  see #18764)
 a slow, exact one that can set the simplex basis (only
InteractiveLPBackend
fits the bill  from #20296)
Ideally, in pure LP mode, both backends would support the basisstatus functions that can transplant the (hopefully) optimal (hopefully)basis from the inexact LP to the exact LP.
If the inexact LP cannot provide a basis (because its "basis" is not a basis due to numerics, or because basisstatus functions are not available), one could at least try to make use of the numerical solution vector and try to reconstruct a basis, like in interiorpointtosimplex crossover (a classical paper: http://www.caam.rice.edu/caam/trs/91/TR9132.pdf)
In MIP mode, could at least try to set the cleanedup numerical solution vector as a known solution, to speed up branchandcut in the exact solver.
Sounds like a big ticket; we'll do this step by step.
#18685 provides the necessary basisstatus functions (for the GLPK backend). #18688 provides a solverindependent interface to these functions. #18804 exposes basis status via backend dictionaries.
Change History (39)
comment:1 Changed 6 years ago by
 Cc dimpase added
comment:2 followup: ↓ 5 Changed 6 years ago by
comment:3 Changed 6 years ago by
 Description modified (diff)
On the other hand, a solverindependent way to get an optimal dual solution is very much welcome, as this is lacking currently, and often needed.
comment:4 Changed 6 years ago by
 Description modified (diff)
comment:5 in reply to: ↑ 2 ; followup: ↓ 6 Changed 6 years ago by
Replying to dimpase:
Is
ppl
(ppl
LP backend, which works with exact arithmetic) too slow for you?
Dima, ppl's implementation of the double description method is very good, but its LP solver is not suitable for problems of even moderate sizes.
comment:6 in reply to: ↑ 5 ; followup: ↓ 7 Changed 6 years ago by
Replying to mkoeppe:
Replying to dimpase:
Is
ppl
(ppl
LP backend, which works with exact arithmetic) too slow for you?Dima, ppl's implementation of the double description method is very good, but its LP solver is not suitable for problems of even moderate sizes.
Would you mind providing an example of PPL choking on an LP doable in exact arithmetic by another solver? We use PPL's LP solver in codesize_upper_bound(...,algorithm="LP")
and never saw a problem... (Although perhaps the difficulty from entry sizes dominate the the one from the dimension in this case).
comment:7 in reply to: ↑ 6 ; followup: ↓ 8 Changed 6 years ago by
Replying to dimpase:
Would you mind providing an example of PPL choking on an LP doable in exact arithmetic by another solver? We use PPL's LP solver in
codesize_upper_bound(...,algorithm="LP")
and never saw a problem... (Although perhaps the difficulty from entry sizes dominate the the one from the dimension in this case).
In our experiments here, we don't actually have numerical difficulties with floatingpoint based solvers; we just want to be sure that we have an exact optimal solution. With #18764 (glp_exact; please review) we have now run some tests to compare performance:
glp_simplex glp_simplex+glp_exact glp_simplex glp_exact +glp_exact ppl + reconstruction in Sage 10 4.20 51.92 7.78 207.07 289.00 11 5.08 58.49 9.43 3451.42 574.72 12 7.55 101.72 11.32 1252.91 808.73 13 7.21 279.08 13.57 1424.28 1019.95 14 8.41 562.97 15.91 1628.54 15 13.10 550.46 18.48 2550.94
As you can see, PPL is much slower than pure glp_exact, and orders of magnitudes slower than glp_simplex followed by glp_exact.
However, currently when we try to reconstruct the solution from the combinatorial basis information, Sage's super slow matrix functions over the rationals get us back to the same order of magnitude as PPL.
It would be interesting to know how the solvers perform on the kind of LPs that you have in mind.
comment:8 in reply to: ↑ 7 Changed 6 years ago by
Replying to mkoeppe:
It would be interesting to know how the solvers perform on the kind of LPs that you have in mind.
LPs I get would be not possible to even enter into a solver without long integers/rationals. That's e.g. behind this function call:
sage: codesize_upper_bound(70,8,2,algorithm="LP") 9695943911863423
more explicitly, you can do
sage: v,p,r=delsarte_bound_hamming_space(70,8,2,return_data=True) sage: p Mixed Integer Program ( maximization, 71 variables, 148 constraints )
constrains of p
have entries as big as 112186277816662845432.
comment:9 Changed 6 years ago by
 Description modified (diff)
 Summary changed from MixedIntegerLinearProgram: Reconstruct exact rational basic solution to MixedIntegerLinearProgram: Reconstruct exact rational/algebraic basic solution
comment:10 Changed 6 years ago by
 Dependencies changed from #18685, #18688 to #18685, #18688, #20296
 Description modified (diff)
 Summary changed from MixedIntegerLinearProgram: Reconstruct exact rational/algebraic basic solution to MixedIntegerLinearProgram/HybridBackend: Reconstruct exact rational/algebraic basic solution
comment:11 Changed 6 years ago by
 Description modified (diff)
comment:12 Changed 6 years ago by
 Branch set to u/mkoeppe/hybrid_backend
comment:13 Changed 5 years ago by
 Commit set to 0b8b78af1c9efbc118513cfde612eccb0bf735a6
 Description modified (diff)
Last 10 new commits:
e2319b5  InteractiveLPBackend.get_variable_value: Guard against standardform transformations

e27f297  InteractiveLPBackend: Make base_ring an init argument

5b0954f  InteractiveLPBackend._variable_type_from_bounds: Add doctests

c4b93aa  InteractiveLPBackend: Fix oldstyle raise statements

b0a3c1c  GenericBackend: Add a missing '# optional  Nonexistent_LP_solver'

3770be0  default_mip_solver: Handle 'InteractiveLP'

d91c776  default_mip_solver, get_solver: Mention InteractiveLP in the documentation

eaede28  get_solver: Add optional base_ring argument

184249d  MixedIntegerLinearProgram: New base_ring init argument

0b8b78a  HybridBackend: first draft

comment:14 Changed 5 years ago by
 Commit changed from 0b8b78af1c9efbc118513cfde612eccb0bf735a6 to 26dad9482713c0a1bc4999b61ce52ed1f109d432
Branch pushed to git repo; I updated commit sha1. This was a forced push. New commits:
26dad94  HybridBackend: first draft

comment:15 Changed 5 years ago by
 Milestone changed from sage6.8 to sage7.4
comment:16 Changed 8 months ago by
 Branch changed from u/mkoeppe/hybrid_backend to u/yzh/hybrid_backend
comment:17 Changed 8 months ago by
 Commit changed from 26dad9482713c0a1bc4999b61ce52ed1f109d432 to 5ee77381d3f8dad0b60137c1f08158ec5668b25e
Branch pushed to git repo; I updated commit sha1. New commits:
5ee7738  change Cython to Python code

comment:18 Changed 8 months ago by
 Commit changed from 5ee77381d3f8dad0b60137c1f08158ec5668b25e to 3431fc483be2286e419a92b6c57e92e701cbffa3
Branch pushed to git repo; I updated commit sha1. New commits:
3431fc4  HybridBackend heritage of GenericBackend

comment:19 Changed 8 months ago by
 Milestone changed from sage7.4 to sage9.3
comment:20 Changed 8 months ago by
 Commit changed from 3431fc483be2286e419a92b6c57e92e701cbffa3 to c1951c9ff301b12d3c5d9fda0bf5ed0a6def282a
Branch pushed to git repo; I updated commit sha1. New commits:
c1951c9  fix docstring style

comment:21 Changed 8 months ago by
 Commit changed from c1951c9ff301b12d3c5d9fda0bf5ed0a6def282a to c15aaac11e8e0650cdf8a1b6f25d4cb5b7b4af37
comment:22 Changed 8 months ago by
 Commit changed from c15aaac11e8e0650cdf8a1b6f25d4cb5b7b4af37 to 44bcb393291d78669ba64ac53f5b490e9ae4de9c
Branch pushed to git repo; I updated commit sha1. This was a forced push. New commits:
3e4e7e9  fix docstrings

69a9775  warmstart interactive_backend.solve() by providing basic_variables

676c74f  warmstart LPProblemStandardForm.run_revised_simplex_method by providing basic_variables

44bcb39  HybridBackend.solve() reconstruct exact solution

comment:23 Changed 8 months ago by
 Status changed from new to needs_review
comment:24 Changed 8 months ago by
comment:25 followup: ↓ 27 Changed 8 months ago by
I don't like the changes to src/sage/numerical/interactive_simplex_method.py
too much.
Instead of extending InteractiveLPProblemStandardForm.run_revised_simplex_method
, I think it's better to use the "run_simplex_method" method of the current_dictionary (a revised dictionary) directly.
Same amount of code, but would avoid making interactive_simplex_method
more complicated.
comment:26 Changed 8 months ago by
 Commit changed from 44bcb393291d78669ba64ac53f5b490e9ae4de9c to 1f9d00adeca30d58c8b23e66570e71d1b7723a4d
comment:27 in reply to: ↑ 25 Changed 8 months ago by
Replying to mkoeppe:
I don't like the changes to
src/sage/numerical/interactive_simplex_method.py
too much.Instead of extending
InteractiveLPProblemStandardForm.run_revised_simplex_method
, I think it's better to use the "run_simplex_method" method of the current_dictionary (a revised dictionary) directly.Same amount of code, but would avoid making
interactive_simplex_method
more complicated.
Good idea. Done.
comment:28 followup: ↓ 33 Changed 8 months ago by
Also, I think it would be better to avoid changing the solve
interface:
 cpdef int solve(self) except 1: + cpdef int solve(self, basic_variables=[]) except 1:
Setting an initial basis could be done in a new method perhaps called warmstart(basic_variables)
or set_basis(basic_variables)
comment:29 followup: ↓ 30 Changed 8 months ago by
And note that generic backend already defines getter methods for the basis status is_variable_basic
, is_variable_nonbasic_at_lower_bound
, is_slack_variable_basic
, is_slack_variable_nonbasic_at_lower_bound
 all of which use variable/row indices.
So this interface:
+ ``basic_variables`` can be one of the following: + +  a list of indices. The indices (starting at 1) correspond to that of the vector formed by `self.interactive_lp_problem().decision_variables()` and `self.interactive_lp_problem().slack_variables()`. Remark that `self.interactive_lp_problem()` can have more variables and constraints than that of `self` if `self` has free variables or `==` constraints. + +  a list of the names of the variables in `self.interactive_lp_problem()`.
looks a bit out of place.
comment:30 in reply to: ↑ 29 Changed 8 months ago by
What would be a better input form? Remove the second option and take only a list of indices?
It is tricky because the indices here correspond to the variables in the interactive_lp_problem().standard_form()
, and they not necessarily the same as the indices for is_variable_basic
etc. which correspond to interactive_lp_problem()
. See the following example.
from sage.numerical.backends.generic_backend import get_solver h = get_solver(solver = ("InteractiveLP")) h.add_variables(1, lower_bound=None, upper_bound=None); h.add_variables(1, lower_bound=0, upper_bound=None); h.add_linear_constraint([(0,2),(1,1)],1,None) h.add_linear_constraint([(0,1),(1,1)],None, 1) h.add_linear_constraint([(0,1),(1,1)],2,2) h.set_objective([0,1]) lp = h.interactive_lp_problem(); view(lp) st = lp.standard_form(); view(st) R = st.coordinate_ring(); R
Replying to mkoeppe:
And note that generic backend already defines getter methods for the basis status
is_variable_basic
,is_variable_nonbasic_at_lower_bound
,is_slack_variable_basic
,is_slack_variable_nonbasic_at_lower_bound
 all of which use variable/row indices.So this interface:
+ ``basic_variables`` can be one of the following: + +  a list of indices. The indices (starting at 1) correspond to that of the vector formed by `self.interactive_lp_problem().decision_variables()` and `self.interactive_lp_problem().slack_variables()`. Remark that `self.interactive_lp_problem()` can have more variables and constraints than that of `self` if `self` has free variables or `==` constraints. + +  a list of the names of the variables in `self.interactive_lp_problem()`.looks a bit out of place.
comment:31 Changed 8 months ago by
 Status changed from needs_review to needs_work
comment:32 Changed 8 months ago by
 Commit changed from 1f9d00adeca30d58c8b23e66570e71d1b7723a4d to dc36f82f2c31e9e0e414289d007b65c0fb4ce700
Branch pushed to git repo; I updated commit sha1. New commits:
dc36f82  set InteractiveLPBackend.current_dictionary and warm start solve from there

comment:33 in reply to: ↑ 28 Changed 8 months ago by
 Status changed from needs_work to needs_review
Replying to mkoeppe:
Also, I think it would be better to avoid changing the
solve
interface: cpdef int solve(self) except 1: + cpdef int solve(self, basic_variables=[]) except 1:Setting an initial basis could be done in a new method perhaps called
warmstart(basic_variables)
orset_basis(basic_variables)
InteractiveLPBackend now has a new method set_dictionary
and a new attribute current_dictionary
. I didn't name it set_basis
, to avoid confusion between the basic variables in a dictionary and those of the generic backend.
comment:34 followup: ↓ 35 Changed 8 months ago by
I still have concerns about the new interface as well. The problem now is that set_basis
is not defined for any backend other than InteractiveLPBackend
(and also its interface refers to internal details of that).
But more importantly, there should be at least one example that illustrates that this new backend does something useful.
comment:35 in reply to: ↑ 34 ; followup: ↓ 36 Changed 8 months ago by
Replying to mkoeppe:
I still have concerns about the new interface as well. The problem now is that
set_basis
is not defined for any backend other thanInteractiveLPBackend
(and also its interface refers to internal details of that).
Do you mean set_dictionary
? It has a different interface than set_basis
(in #18688 to be written) since the interactive_simplex_method
is indeed very special in running the simplex method on the standard form where more auxiliary variables are introduced. I couldn't think of a way to unify InteractiveLPBackend.set_dictionary
with other backends' set basis methods.
But more importantly, there should be at least one example that illustrates that this new backend does something useful.
True. The current implementation enables what the first paragraph of this ticket states, by using solver = ("GLPK", "InteractiveLP")
as in the examples in solve
. However, it is hard to argue with doctests that this new hybrid backend is useful. Shall we close the ticket with the tag "wontfix" as in #18804?
comment:36 in reply to: ↑ 35 Changed 8 months ago by
Replying to yzh:
Replying to mkoeppe:
I still have concerns about the new interface as well. The problem now is that
set_basis
is not defined for any backend other thanInteractiveLPBackend
(and also its interface refers to internal details of that).Do you mean
set_dictionary
?
Yes, that's what I meant.
comment:37 Changed 7 months ago by
 Milestone changed from sage9.3 to sage9.4
Setting new milestone based on a cursory review of ticket status, priority, and last modification date.
comment:38 Changed 6 months ago by
 Commit changed from dc36f82f2c31e9e0e414289d007b65c0fb4ce700 to 71ae94d57b06bd58a1cf02a9b286bcbe881a0fcc
Branch pushed to git repo; I updated commit sha1. New commits:
71ae94d  Replace not hasattr(self, 'xyz') by self.xyz is None

comment:39 Changed 2 months ago by
 Milestone changed from sage9.4 to sage9.5
Setting a new milestone for this ticket based on a cursory review.
Is
ppl
(ppl
LP backend, which works with exact arithmetic) too slow for you?