Open TheCoolRob opened 6 years ago
I would like to know more. As we are looking at MIP cuts and presolve in the context of the CP-SAT solver, this is an interesting topic.
Which cuts are you adding? Laurent Perron | Operations Research | lperron@google.com | (33) 1 42 68 53 00
Le mar. 30 oct. 2018 à 09:31, TheCoolRob notifications@github.com a écrit :
Hi
Following on from #666 https://github.com/google/or-tools/issues/666 - I found the reason why I wasn't using the CPModelProto. I've got some custom cutting planes I'm implementing over a SAT model and using the CutGenerator callbacks.
I wasn't able to find an interface on the CPModelProto that allows for adding a cut-callback. I would prefer to work at the CPModel level rather than SAT (as per Laurents' suggestion) but perhaps I'm just missing the part of the API that gives the appropriate functionality?
It could be that it's not available at that level - passing a pointer to a callback over the proto model seems like a bad idea and something that probably wouldn't be supported. I suppose this issue is simply then for confirmation that this is the case or not.
Thanks
— You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub https://github.com/google/or-tools/issues/912, or mute the thread https://github.com/notifications/unsubscribe-auth/AKj17UMUx_8daiyrDs_wViLDfU7RHi5jks5uqA5ogaJpZM4YBaXJ .
Hi Laurent. So I was using the SAT cut generator layer that was available and extended your circuit cuts that were written.
Nice!
Do you use the two commodity flow as a relaxation of the TSP?
The current flow with the cp_model API is that the constraint exists in the proto. When parsing and loading the model, we extract it in the CP-SAT layer, and in the LP.
Thus, improvements to the circuit falls directly into this schema.
Are comb, blossom, connect cuts related to the circuit constraint too? Do you have a reference for them?
Thanks
Laurent Perron | Operations Research | lperron@google.com | (33) 1 42 68 53 00
Le mar. 6 nov. 2018 à 09:27, TheCoolRob notifications@github.com a écrit :
Hi Laurent. So I was using the SAT cut generator layer that was available and extended your circuit cuts that were written.
- The SAT cut gen only permitted cutting the root node once. It was relatively straight forward to modify your subtour elimination constraints (sec's) to run recursively until the min-cut returns cuts that already exist. This typically creates a tighter model for circuits - I got a nice improvement from that but forgot to push the change on the public repo (I wasn't sure how it would perform in all situations).
- I added comb, blossom and connect cuts - but because of the representation issues in SAT (complete graph, directed edges) I moved to a two-commodity flow representation in MIP which gave me pretty tight bounds for the problem size I was tackling (mapping the undirected edge formulation cuts into the commodity flow representation).
- I wanted to scale it up and use SAT for local search under the MIP but realised that without a cut-callback on the CPModelProto I wouldn't be able to generate cuts with a common definition of the cuts I'm using at the moment. It produced an itch that resulted in this github issue.. Thanks
— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/google/or-tools/issues/912#issuecomment-436169536, or mute the thread https://github.com/notifications/unsubscribe-auth/AKj17XFDC0Uol0UrS4UyVuxFWxeVGl5pks5usUfvgaJpZM4YBaXJ .
This is a good reference for the cuts I've described related to the TSP: https://www.math.uwaterloo.ca/~bico/papers/dfj_mathprog.pdf
I'm using a relaxed version of the two-commodity flow, but with column generation for the edges (my comment about not needing a complete graph, makes the models much more compact). There several other cuts that are specific to the problem outside of the circuit feature (like capacity cuts and precedence cuts) as it's effectively a pickup and dropoff problem on top of the standard TSP. You can have both models side by side (flow variables which induce binary variables in the circuit representation). This gives you the best of both worlds, the tight bound from the DFJ representation and cumul dimension cuts from the flow model. Requires a bit of effort to keep them in sync but it's a good approach - as long as you have a cut callback at every node ;-)
Since we're on the topic of model representations - I saw you made a comment on the forums regarding upcoming changes to the routing library (re it's deprecation). I've seen some hints in the code that point to a more precise approach that you may be planning. How drastic are the changes you're planning and do you have provisional timeline for the release?
For routing update, once we migrate our github source code to use the abseil-cpp
third party dependency, then we will also be able to get in sync with our internal source code ~1 year ahead concerning the routing library.
You can follow advancement at #871 or looking at the Work In Progress abseil
branch directly.
note: commit contains both abseil-cpp migration and backport of the "new" routing library
Hi,
thanks for the pointers. So you are trying to solve a PDP with capacity (and maybe time windows) exactly using SAT. I like that.
The routing library is not under maintenance, only the original CP code will only be extended to support the routing part. I strongly recommend persons not using the routing library to use the CP-SAT solver. The new API is visible in the abseil branch of or-tools on github:
Here is the ported example:
https://github.com/google/or-tools/blob/abseil/examples/cpp/tsp.cc
3 notable changes:
I ported all examples to the new API in very little time, and the benefits in term of code (no more mixing node and indices) are really worth it. Laurent Perron | Operations Research | lperron@google.com | (33) 1 42 68 53 00
Le mar. 6 nov. 2018 à 10:56, TheCoolRob notifications@github.com a écrit :
This is a good reference for the cuts I've described related to the TSP: https://www.math.uwaterloo.ca/~bico/papers/dfj_mathprog.pdf
I'm using a relaxed version of the two-commodity flow, but with column generated for the edges (my comment about not needing a complete graph, makes the models much more compact). There several other cuts that are specific to the problem outside of the circuit feature (like capacity cuts and precedence cuts) as it's effectively a pickup and dropoff problem on top of the standard TSP. You can have both models side by side (flow variables which induce binary variables in the circuit representation). This gives you the best of both worlds, the tight bound from the DFJ representation and cumul dimension cuts from the flow model. Requires a bit of effort to keep them in sync but it's a good approach - as long as you have a cut callback at every node ;-)
Since we're on the topic of model representations - I saw you made a comment on the forums regarding upcoming changes to the routing library (re it's deprecation). I've seen some hints in the code that point to a more precise approach that you may be planning. How drastic are the changes you're planning and do you have provisional timeline for the release?
— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/google/or-tools/issues/912#issuecomment-436195225, or mute the thread https://github.com/notifications/unsubscribe-auth/AKj17audFbuQnjuT-jFPVGyFjxlf7ByQks5usVzPgaJpZM4YBaXJ .
Can you share the code?
We will reimplement it, and as David works at Google, we can ask questions directly, but it will guide us.
Thanks Laurent Perron | Operations Research | lperron@google.com | (33) 1 42 68 53 00
Le mar. 6 nov. 2018 à 11:08, Laurent Perron lperron@google.com a écrit :
Hi,
thanks for the pointers. So you are trying to solve a PDP with capacity (and maybe time windows) exactly using SAT. I like that.
The routing library is not under maintenance, only the original CP code will only be extended to support the routing part. I strongly recommend persons not using the routing library to use the CP-SAT solver. The new API is visible in the abseil branch of or-tools on github:
Here is the ported example:
https://github.com/google/or-tools/blob/abseil/examples/cpp/tsp.cc
3 notable changes:
- we do create an index manager that sync nodes in the routing graph and indices in the underlying model.
- distance and capacity callbacks are registered on the routing model, and we use their index (returned by the register method) in dimension, costs... instead of the callback.
- in non C++ languages, callbacks inherits from IntIntToLong instead of NodeEvaluator2
I ported all examples to the new API in very little time, and the benefits in term of code (no more mixing node and indices) are really worth it. Laurent Perron | Operations Research | lperron@google.com | (33) 1 42 68 53 00
Le mar. 6 nov. 2018 à 10:56, TheCoolRob notifications@github.com a écrit :
This is a good reference for the cuts I've described related to the TSP: https://www.math.uwaterloo.ca/~bico/papers/dfj_mathprog.pdf
I'm using a relaxed version of the two-commodity flow, but with column generated for the edges (my comment about not needing a complete graph, makes the models much more compact). There several other cuts that are specific to the problem outside of the circuit feature (like capacity cuts and precedence cuts) as it's effectively a pickup and dropoff problem on top of the standard TSP. You can have both models side by side (flow variables which induce binary variables in the circuit representation). This gives you the best of both worlds, the tight bound from the DFJ representation and cumul dimension cuts from the flow model. Requires a bit of effort to keep them in sync but it's a good approach - as long as you have a cut callback at every node ;-)
Since we're on the topic of model representations - I saw you made a comment on the forums regarding upcoming changes to the routing library (re it's deprecation). I've seen some hints in the code that point to a more precise approach that you may be planning. How drastic are the changes you're planning and do you have provisional timeline for the release?
— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/google/or-tools/issues/912#issuecomment-436195225, or mute the thread https://github.com/notifications/unsubscribe-auth/AKj17audFbuQnjuT-jFPVGyFjxlf7ByQks5usVzPgaJpZM4YBaXJ .
@Mizux thanks for the heads up on the branch - does that mean that branch is currently up to date with the internal state of the routing solver or that it be updated once the abseil move is complete?
@lperron thanks for the link to the sample. Moving away from those two layers of indexes will mean a 50% reduction in questions on the forum :1st_place_medal: . I suppose my question was more about the structure of the way in which the solver works - specifically, pickup and dropoff models (from my testing at any rate) don't seem to produce good solutions. Some of the seeding mechanisms just aren't designed to handle the dependencies between orders and the local search heuristics chain together in a really inefficient way (mostly O(n^2) local search operators for pairs). This isn't a fair complaint, because the routing layer does a lot really well, but it felt like perhaps if you were tackling some restructuring of the design then I was hoping this would be on the list?
Lastly on the code - I can't share my implementation unfortunately (its a commercial implementation) - but luck is on your side because you can ask David if you can port the concorde procedures (the canonical implementation referenced in that paper) into the open source world (currently it's only available under for academic usage). If he would prefer that you re-write them I'd be happy to help contribute in this regard if/when you get around to it as there is no publicly available unrestricted-use license for those cut procedures (not to my knowledge, happy to be corrected if wrong).
@TheCoolRob Should be in sync modulo 2 weeks I would say....
A lot of work has been done on the routing library to improve performance of the local search, and the first solution heuristics w.r.t. extra constraints like PDP.
I would expect this version to be much better.
About code sharing, I understand. If there are anything you can share, please send it to me.
Thanks Laurent Perron | Operations Research | lperron@google.com | (33) 1 42 68 53 00
Le mar. 6 nov. 2018 à 12:45, Mizux notifications@github.com a écrit :
@TheCoolRob https://github.com/TheCoolRob Should be in sync modulo 2 weeks I would say....
— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub https://github.com/google/or-tools/issues/912#issuecomment-436225165, or mute the thread https://github.com/notifications/unsubscribe-auth/AKj17bGwMYKEI1vh1JGsvhSrrJI2MkyLks5usXY9gaJpZM4YBaXJ .
Thanks @lperron. Looking forward to testing the new release when it's available with these changes.
Looping back to the original topic - a cut-call back that can be specified through the CPModelProto layer - does this idea make sense in your design for the CP-SAT layer or is this a different kind of animal?
Interesting idea. We need to think of the architecture. In the meantime, you could build from source and modify the source locally. If you have code you can share, we can integrate it/reimplement it, or add an API for callbacks.
--Laurent Laurent Perron | Operations Research | lperron@google.com | (33) 1 42 68 53 00
Le mar. 6 nov. 2018 à 14:57, TheCoolRob notifications@github.com a écrit :
Thanks @lperron https://github.com/lperron. Looking forward to testing the new release when it's available with these changes.
Looping back to the original topic - a cut-call back that can be specified through the CPModelProto layer - does this idea make sense in your design for the CP-SAT layer or is this a different kind of animal?
— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub https://github.com/google/or-tools/issues/912#issuecomment-436261097, or mute the thread https://github.com/notifications/unsubscribe-auth/AKj17bVwwiPmXarVyYNVh0SCOrogULzUks5usZVWgaJpZM4YBaXJ .
I suspect callback for the api would be quite simple since a cut is only a function of an incumbent. So a solution vector in, vector of inequalities out. That way you're also decoupled from my situation of having 2 different representations of the same model (i.e. the user is responsible for the translation).
I agree that some working code is a sensible way forward. I'll put together a use-case and perhaps it will be clearer how to extract it out into something that can be used generically.
Yes and no, we have presolve in between. This can become tricky. Laurent Perron | Operations Research | lperron@google.com | (33) 1 42 68 53 00
Le mar. 6 nov. 2018 à 15:12, TheCoolRob notifications@github.com a écrit :
I suspect callback for the api would be quite simple since a cut is only a function of an incumbent. So a solution vector in, vector of inequalities out. That way you're also decoupled from my situation of having 2 different representations of the same model (i.e. the user is responsible for the translation).
I agree that some working code is a sensible way forward. I'll put together a use-case and perhaps it will be clearer how to extract it out into something that can be used generically.
— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub https://github.com/google/or-tools/issues/912#issuecomment-436265884, or mute the thread https://github.com/notifications/unsubscribe-auth/AKj17YfAasaewBh5eRhc5hSutI7YqOg2ks5usZi2gaJpZM4YBaXJ .
Presolve does complicate it - I forgot about that. It feels like there's no simple way around presolve without disabling it because a cut might then invalidate a presolve aggregation that would need be unrolled? Tricky might be an understatement - I can think of examples where many variables might need to be added and constraints adjusted (removed and relinked) in a model. I'll pretend I didn't see this comment then and just disable presolve :joy:
What we do is presolve produces a new model proto, plus a channeling model proto. We solve the inner model, and populate the solution of the initial model with the solution of the internal model, and the channeling model.
Thus, the right framework is to take a write a cut generator for the constraint in the model. This being said, it is not impossible to add a constraint that is only build by presolve/aggregation, and thus extracted to the LP with more information. But we try to keep it simple.
Thanks Laurent Perron | Operations Research | lperron@google.com | (33) 1 42 68 53 00
Le mar. 6 nov. 2018 à 16:16, TheCoolRob notifications@github.com a écrit :
Presolve does complicate it - I forgot about that. It feels like there's no simple way around presolve without disabling it because a cut might then invalidate a presolve aggregation that would need be unrolled? Tricky might be an understatement - I can think of examples where many variables might need to be added and constraints adjusted (removed and relinked) in a model. I'll pretend I didn't see this comment then and just disable presolve 😂
— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub https://github.com/google/or-tools/issues/912#issuecomment-436288566, or mute the thread https://github.com/notifications/unsubscribe-auth/AKj17Q76u-xctgaKTbladncMroyIEVl5ks5usafPgaJpZM4YBaXJ .
Hi Laurent,
This is a follow up question to this implementation. In a standard TSP formulation one will only work with undirected edges. I see that the SAT-circuit formulation seems to have directed arcs (or edges, I think it's largely a discipline preferences). I think it's a lot easier to work with directed edges because it means one's intuition of how to walk the graph is modelled exactly the same way. That said, there are quite a few drawbacks when it comes to graph operations (such as chained Lin-Kernighan) which can't be performed in a computationally efficient way on a directed graph (edge flips then require propagation through the whole circuit to validate rather than on the local exchange - so O(n) for each operation rather than approximately O(1)).
There is a standard Asymmetric -> Symmetric TSP transformation which just adds n fixed undirected edges to the graph and allows one to treat the asymmetric problem as a symmetric one (getting the speed up from the independence of the local exchanges). It requires that you can handle a sparse set of edges (which I think you already support in the circuit constraint). It's a bit complicated to map everything into the circuit layer from an undirected context into a directed context - so I suppose I'm asking if you've got any plans in this regard or whether I should forge on?
Right now, we are lagging behind you. We currently have no plan to lift ATSP to TSP. Therefore, you can forge ahead.
One thing I do not know is whether the symmetric TSP cuts are valid w.r.t. additional constraints, or if they are dominance rules valid only on the pure problem. I need a circuit constraint in the TSP to solve complex routing (vrp, pdptw, routing mixed with scheduling or assignment...) problems. I do not know if the comb, blossom cuts are useful in that setting.
Thanks
Laurent Perron | Operations Research | lperron@google.com | (33) 1 42 68 53 00
Le lun. 12 nov. 2018 à 08:34, TheCoolRob notifications@github.com a écrit :
Hi Laurent,
This is a follow up question to this implementation. In a standard TSP formulation one will only work with undirected edges. I see that the SAT-circuit formulation seems to have directed arcs (or edges, I think it's largely a discipline preferences). I think it's a lot easier to work with directed edges because it means one's intuition of how to walk the graph is modelled exactly the same way. That said, there are quite a few drawbacks when it comes to graph operations (such as chained Lin-Kernighan) which can't be performed in a computationally efficient way on a directed graph (edge flips then require propagation through the whole circuit to validate rather than on the local exchange - so O(n) for each operation rather than approximately O(1)).
There is a standard Asymmetric -> Symmetric TSP transformation which just adds n fixed undirected edges to the graph and allows one to treat the asymmetric problem as a symmetric one (getting the speed up from the independence of the local exchanges). It requires that you can handle a sparse set of edges (which I think you already support in the circuit constraint). It's a bit complicated to map everything into the circuit layer from an undirected context into a directed context - so I suppose I'm asking if you've got any plans in this regard or whether I should forge on?
— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub https://github.com/google/or-tools/issues/912#issuecomment-437783641, or mute the thread https://github.com/notifications/unsubscribe-auth/AKj17YrQBQJpzcDMpZJ1RJ56j_YZ78Ucks5uuSSPgaJpZM4YBaXJ .
I will forge ahead then. The mapping is a bit painful because you have to re-map all the equations and adjust coefficients - I'm almost certainly going to lose a day to a bug somewhere :joy:
In terms of the validity of the cuts. They are still valid in the circuit context and they exist in addition to the other constraints you've mentioned here. In the implementation I mentioned earlier, I'm doing complex pickup/dropoff with all the standard TSP cuts but also Precedence Cuts, Capacity Cuts, etc. The TSP cuts do 95% of the work, the side constraints provide cuts that close the gap to around 1% in most cases. Without all the cuts you typically have relaxations that approximate the minimum spanning tree (which leaves a gap of around 50%). In my opinion, it's worth the effort, but my use case is very specific to routing.
There are some assumptions in the tsp cuts (non-negativity of edge weights) and I think the graph needs to be metric (adhere to the triangle inequality) - which can be tricky if the routing algorithms aren't using admissible heuristics on a road network. This is somewhat of a side point though, because one can detect if these requirements apriori and disable the affected cuts if the assumptions aren't met.
Sorry this response is so overdue - needed to carve out a morning to write this up. I've wired up the exact search directly in the SAT model (so bypassing the CPModelProto
layer) and found a few interesting things (which you probably already knew).
The test problem is a simple TSP with 100 nodes on a asymmetric distance matrix. I previously tested this with CBC as the backing solver and am able to confirm an optimal in around 13.5 seconds. The MIP model uses the standard Asymmetric->Symmetric transformation. The procedure is a little hacky (because CBC doesn't support non-root node cuts) so the output below is just the last solution.
Clp0032I Optimal objective 1527.821254 - 504 iterations time 0.042
0 solution(s) were found (by branching) at an average depth of 0
average value of variable being branched on was 0.493394
101 nodes were cutoff at an average depth of 9.80198 with iteration count of 34.6931
100 nodes were solved at an average depth of 8.75 with iteration count of 24.22
Down 101 nodes (66 first, 35 second) - 61 cutoff, rest decrease numinf 1.65347 increase obj 0.0258301
Up 100 nodes (35 first, 65 second) - 40 cutoff, rest decrease numinf 2.76 increase obj 0.243179
Total solve time for 100 nodes: 13.5156
It takes 9 relaxation iterations and a bit of B&B to get to that solution.
I set up the same problem against the SAT layer and got the following (using a directed formulation with the LP relaxation cut generator):
I1210 21:53:50.901527 24446 optimization.cc:1131] #Boolean_variables:10000
status: OPTIMAL
objective: 1532.306
best_bound: 1532.306
booleans: 10000
conflicts: 3974
branches: 26475
propagations: 652130
walltime: 16.362212
deterministic_time: 0.343003
SAT solve time: 0
Total solve time for 100 points: 19.5197
The cut generation time is around 3 seconds of the total time in the SAT solve.
I scaled the objective to 3 decimal places and confirmed the solution is identical in both. The optimal solution reported in the cbc solution is a bit lower - when you project is back up into a directed formulation it collapses to the same solution as produced by the SAT-solver. So the SAT apporach is a bit slower than a more traditional MIP approach (although the CBC branching strategy is quite sophisticated - if I use GLPK (which is quicker - but less smart) it takes around twice as long to complete the optimality proof).
I'm using the MinimizeIntegerVariableWithLinearScanAndLazyEncoding
+ ExploitIntegerLpSolution
+ UnassignedVarWithLowestMinAtItsMinHeuristic
approach in the SAT solver. I'm not sure if you have any recommendations on a potentially better configuration?
I ran a series of small problems use the CPModelProto
layer with the default strategy - this is what the growth of computation time looked like:
It's not surprising that the growth is so aggressive (giving it is a factorial model) - but I suppose what is quite nice is that it's clear what the benefit of the cuts are in the model. You're able to go from a model of size +-22 to a size of 100 (which is 10^136 larger in search space) with an optimality proof :-). I plotted the performance of the SAT model with the cuts on the red line.
So I suppose the original question on this github issue was whether it would be possible to have some kind of call-back procedure to add cuts - but having just been through this implementation I almost feel like it's a somewhat unreasonable request. The only solutions I could come up with are very cpp specific and probably won't work in other lanugages (perhaps this is an advanced enough feature that this would be okay? - untill someone else demands it I suppose).
I'm going to do a bit of digging to understand how to get the volatility down on the SAT solver with the cut generator (the red line). It would seem there are some easy-win paths which it follows occassionally which would be good to heuristically pick up on if possible. The MIP solver doesn't have this kind of volatility in it's branching strategy - so perhaps it makes sense to use something other than UnassignedVarWithLowestMinAtItsMinHeuristic
.
Overall I'm satisfied - I just felt like I owed a bit of upstream feedback. Hope it's useful.
Thanks,
how do you model the TSP with the CpModelProto interface? Do you use the circuit constraint (which would enable cuts) ?
In our test, we scale aTSP to ~100-150 nodes.
Laurent Perron | Operations Research | lperron@google.com | (33) 1 42 68 53 00
Le mar. 11 déc. 2018 à 08:07, TheCoolRob notifications@github.com a écrit :
Sorry this response is so overdue - needed to carve out a morning to write this up. I've wired up the exact search directly in the SAT model (so bypassing the CPModelProto layer) and found a few interesting things (which you probably already knew).
The test problem is a simple TSP with 100 nodes on a asymmetric distance matrix. I previously tested this with CBC as the backing solver and am able to confirm an optimal in around 13.5 seconds. The MIP model uses the standard Asymmetric->Symmetric transformation. The procedure is a little hacky (because CBC doesn't support non-root node cuts) so the output below is just the last solution.
Clp0032I Optimal objective 1527.821254 - 504 iterations time 0.042 0 solution(s) were found (by branching) at an average depth of 0 average value of variable being branched on was 0.493394 101 nodes were cutoff at an average depth of 9.80198 with iteration count of 34.6931 100 nodes were solved at an average depth of 8.75 with iteration count of 24.22 Down 101 nodes (66 first, 35 second) - 61 cutoff, rest decrease numinf 1.65347 increase obj 0.0258301 Up 100 nodes (35 first, 65 second) - 40 cutoff, rest decrease numinf 2.76 increase obj 0.243179 Total solve time for 100 nodes: 13.5156
It takes 9 relaxation iterations and a bit of B&B to get to that solution.
I set up the same problem against the SAT layer and got the following (using a directed formulation with the LP relaxation cut generator):
I1210 21:53:50.901527 24446 optimization.cc:1131] #Boolean_variables:10000 status: OPTIMAL objective: 1532.306 best_bound: 1532.306 booleans: 10000 conflicts: 3974 branches: 26475 propagations: 652130 walltime: 16.362212 deterministic_time: 0.343003 SAT solve time: 0 Total solve time for 100 points: 19.5197
The cut generation time is around 3 seconds of the total time in the SAT solve.
I scaled the objective to 3 decimal places and confirmed the solution is identical in both. The optimal solution reported in the cbc solution is a bit lower - when you project is back up into a directed formulation it collapses to the same solution as produced by the SAT-solver. So the SAT apporach is a bit slower than a more traditional MIP approach (although the CBC branching strategy is quite sophisticated - if I use GLPK (which is quicker - but less smart) it takes around twice as long to complete the optimality proof).
I'm using the MinimizeIntegerVariableWithLinearScanAndLazyEncoding + ExploitIntegerLpSolution + UnassignedVarWithLowestMinAtItsMinHeuristic approach in the SAT solver. I'm not sure if you have any recommendations on a potentially better configuration?
I ran a series of small problems use the CPModelProto layer with the default strategy - this is what the growth of computation time looked like: [image: image] https://user-images.githubusercontent.com/3518425/49783799-ba88a480-fd23-11e8-82a5-1f8df9b7d595.png
It's not surprising that the growth is so aggressive (giving it is a factorial model) - but I suppose what is quite nice is that it's clear what the benefit of the cuts are in the model. You're able to go from a model of size +-22 to a size of 100 (which is 10^136 larger in search space) with an optimality proof :-). I plotted the performance of the SAT model with the cuts on the red line.
So I suppose the original question on this github issue was whether it would be possible to have some kind of call-back procedure to add cuts - but having just been through this implementation I almost feel like it's a somewhat unreasonable request. The only solutions I could come up with are very cpp specific and probably won't work in other lanugages (perhaps this is an advanced enough feature that this would be okay? - untill someone else demands it I suppose).
I'm going to do a bit of digging to understand how to get the volatility down on the SAT solver with the cut generator (the red line). It would seem there are some easy-win paths which it follows occassionally which would be good to heuristically pick up on if possible. The MIP solver doesn't have this kind of volatility in it's branching strategy - so perhaps it makes sense to use something other than UnassignedVarWithLowestMinAtItsMinHeuristic.
Overall I'm satisfied - I just felt like I owed a bit of upstream feedback. Hope it's useful.
— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub https://github.com/google/or-tools/issues/912#issuecomment-446095392, or mute the thread https://github.com/notifications/unsubscribe-auth/AKj17Sy9e3TuaJhHVznScoV4jfqLNtzjks5u31nEgaJpZM4YBaXJ .
This is a snippet of the CPModelProto that I'm using:
CpModelProto cp_model;
auto new_variable = [&cp_model](int64 lb, int64 ub) {
CHECK_LE(lb, ub);
const int index = cp_model.variables_size();
IntegerVariableProto *const var = cp_model.add_variables();
var->add_domain(lb);
var->add_domain(ub);
return index;
};
auto new_linear_constraint =
[&cp_model](vector<int> indexes, vector<int> coef, double lb, double ub) {
LinearConstraintProto *const lin =
cp_model.add_constraints()->mutable_linear();
for (int i = 0; i < indexes.size(); i++) {
lin->add_vars(indexes[i]);
lin->add_coeffs(coef[i]);
}
lin->add_domain(lb);
lin->add_domain(ub);
};
vector<vector<int>> vars;
auto objective = cp_model.mutable_objective();
for (int i = 0; i < basedim; i++) {
vector<int> tmp;
for (int j = 0; j < basedim; j++) {
int ub = i == j ? 0 : 1; // diagonals not allowed
tmp.push_back(new_variable(0, ub));
objective->add_vars(tmp.back());
objective->add_coeffs((int64)(distance_matrix[i][j] * 1000));
}
vars.push_back(tmp);
}
for (int i = 0; i < basedim; i++) {
vector<int> indexes;
vector<int> coefficients;
for (int j = 0; j < basedim; j++) {
indexes.push_back(vars[i][j]);
coefficients.push_back(1.0);
}
new_linear_constraint(indexes, coefficients, 1, 1);
}
auto circuit = cp_model.add_constraints()->mutable_circuit();
for (int i = 0; i < basedim; i++) {
vector<int> indexes;
vector<int> coefficients;
for (int j = 0; j < basedim; j++) {
indexes.push_back(vars[j][i]);
coefficients.push_back(1.0);
circuit->add_tails(i);
circuit->add_heads(j);
circuit->add_literals(vars[i][j]);
}
new_linear_constraint(indexes, coefficients, 1, 1);
}
Model model;
const CpSolverResponse response = SolveCpModel(cp_model, &model);
std::cout << "SAT solve time: " << response.wall_time() << std::endl;
This doesn't use the connect cuts that you've written up for the model from what I can tell (unless there's some magic behind the circuit when the model gets handled behind the proto layer).
I can run a separate test against your connect cuts to confirm where it lies on the performance curve?
OK, you need to add linearization_level:2 to the sat parameters.
SatParameters parameters; parameters.set_linearization_level(2); model.Add(NewSatParameters(parameters));
Laurent Perron | Operations Research | lperron@google.com | (33) 1 42 68 53 00
Le mar. 11 déc. 2018 à 08:49, TheCoolRob notifications@github.com a écrit :
This is a snippet of the CPModelProto that I'm using:
CpModelProto cp_model;
auto new_variable = [&cp_model](int64 lb, int64 ub) { CHECK_LE(lb, ub); const int index = cp_model.variables_size(); IntegerVariableProto *const var = cp_model.add_variables(); var->add_domain(lb); var->add_domain(ub); return index; };
auto new_linear_constraint = [&cp_model](vector
indexes, vector coef, double lb, double ub) { LinearConstraintProto *const lin = cp_model.add_constraints()->mutable_linear(); for (int i = 0; i < indexes.size(); i++) { lin->add_vars(indexes[i]); lin->add_coeffs(coef[i]); } lin->add_domain(lb); lin->add_domain(ub); }; vector<vector
> vars; auto objective = cp_model.mutable_objective();
for (int i = 0; i < basedim; i++) { vector
tmp; for (int j = 0; j < basedim; j++) { int ub = i == j ? 0 : 1; // diagonals not allowed tmp.push_back(new_variable(0, ub)); objective->add_vars(tmp.back()); objective->add_coeffs((int64)(distance_matrix[i][j] * 1000)); } vars.push_back(tmp); } for (int i = 0; i < basedim; i++) { vector
indexes; vector coefficients; for (int j = 0; j < basedim; j++) { indexes.push_back(vars[i][j]); coefficients.push_back(1.0); } new_linear_constraint(indexes, coefficients, 1, 1); } auto circuit = cp_model.add_constraints()->mutable_circuit();
for (int i = 0; i < basedim; i++) { vector
indexes; vector coefficients; for (int j = 0; j < basedim; j++) { indexes.push_back(vars[j][i]); coefficients.push_back(1.0); circuit->add_tails(i); circuit->add_heads(j); circuit->add_literals(vars[i][j]); } new_linear_constraint(indexes, coefficients, 1, 1);
}
Model model; const CpSolverResponse response = SolveCpModel(cp_model, &model); std::cout << "SAT solve time: " << response.wall_time() << std::endl;
This doesn't use the connect cuts that you've written up for the model from what I can tell (unless there's some magic behind the circuit when the model gets handled behind the proto layer).
I can run a separate test against your connect cuts to confirm where it lies on the performance curve?
— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub https://github.com/google/or-tools/issues/912#issuecomment-446104228, or mute the thread https://github.com/notifications/unsubscribe-auth/AKj17Q-3lbtskTIJnIAqff_s0WG0v0qYks5u32OBgaJpZM4YBaXJ .
That makes a big difference.
I think you've probably optimised quite of a few of the operations between the LP and SAT search procedures in the CpModelProto layer because the search times seem a bit more predictable with this search procedure.
I think one way for me to use this mode of running is just to preprocess the cuts I suspect the model will require and provide them up front. I'm going to run a few tests and see if this helps speed up the search in general or not (I suspect there is a tiny overhead for each cut which a fast SAT solve might implicitly be working around).
The other alternative is to add cut methods directly where your other cuts are in linear_programming_constraint.cc
and test that way around.
So I ran the setting you recommended against my data set to get a comparison. The volatility is definitely lower but it caps out at around 90 points (I left it running for about 3 hours before killing the process) - the models were taking progressively longer to solve.
I created a model which uses a mixture of the CpModelProto
and the Model
class. So pre-populating the variables in the model and setting up the LP class linked to the Model - then calling
const CpSolverResponse response = SolveCpModel(cp_model, &model);
with the linearisation settings enabled as you suggested. This is an attempt to blend the cut generators and the CPProtoModel solve performance with your linearisation settings. It's calling back to the cut generator and the cuts generate are being passed up to the LP correctly. My concern though is around performance - for some reason, trying it this way is much slower than passing the whole CpModelProto down in a single chunk. I'm worried there are either poorly linked variable relations (potentially bypassing preprocessing steps) or duplicate variables by doing it this way. I did some basic testing before writing up the whole algorithm and the moment I added integer variables between the two models there was a 400x performance hit (boolean vars were fine).
It seems like a better approach would be to bypass the CpModelProto
layer and rather replicate the outline of the underlying default search for this problem type so that I can do an apples with apples comparison, and having typed this, I realise that was the suggestion to myself in my previous comment. I'll give that a try - unless there is a simpler alternative?
A quick update on this. I added a custom callback to the Model
object and implemented a cut generator callback in the TryToAddCutGenerators
like this:
auto cb = m->GetOrCreate<CutGenCallback>();
if (cb->IsValid()) {
cut_generators->push_back(
cb->GetCallbackGen(num_nodes, tails, heads, vars));
} else {
cut_generators->push_back(CreateStronglyConnectedGraphCutGenerator(
num_nodes, tails, heads, vars));
}
It's a little dirty but it did the job. I had to do a bit of realignment of the model after presolve made some reasonable changes but nothing too complex luckily. I ran a test on the problem I'm working with (which has some very tricky choice points) - as I mentioned in my previous comment, sometimes the Default SAT search goes on a bit of a tangent, so I put in a time limit of 200 seconds. It still makes for an interesting contrast for what the additional TSP cuts do for the model (red line) over the default SAT search with the recommended settings (black line):
I used a default smoother to show the approximate trends. I suspect that this approach could be improved on further (I only implemented 2 addition cut types) - also the LP exploitation logic or variable selection would probably help reduce the volatility I'm seeing in the runtimes (just to get the search a bit more predictable) although that's probably a bit less general and more circuit specific.
I'm not sure if you've discussed any thoughts on this internally?
I also noticed you've started working on an exact search for the routing problem (there's a CVRP cut generator down there) - have you run any performance comparisons between this and the best known approaches?
Thanks a lot.
We are following this very closely. We are thinking about improving out code and our cuts.
It just takes time. But this helps tremendously.
Thanks again. Laurent Perron | Operations Research | lperron@google.com | (33) 1 42 68 53 00
Le ven. 14 déc. 2018 à 16:35, TheCoolRob notifications@github.com a écrit :
A quick update on this. I added a custom callback to the Model object and implemented a cut generator callback in the TryToAddCutGenerators like this:
auto cb = m->GetOrCreate<CutGenCallback>(); if (cb->IsValid()) { cut_generators->push_back( cb->GetCallbackGen(num_nodes, tails, heads, vars)); } else { cut_generators->push_back(CreateStronglyConnectedGraphCutGenerator( num_nodes, tails, heads, vars)); }
It's a little dirty but it did the job. I had to do a bit of realignment of the model after presolve made some reasonable changes but nothing too complex luckily. I ran a test on the problem I'm working with (which has some very tricky choice points) - as I mentioned in my previous comment, sometimes the Default SAT search goes on a bit of a tangent, so I put in a time limit of 200 seconds. It still makes for an interesting contrast for what the additional TSP cuts do for the model (red line) over the default SAT search with the recommended settings (black line): [image: image] https://user-images.githubusercontent.com/3518425/50011948-05c1e200-ffc6-11e8-85ba-a767af4cd14e.png
I used a default smoother to show the approximate trends. I suspect that this approach could be improved on further (I only implemented 2 addition cut types) - also the LP exploitation logic or variable selection would probably help reduce the volatility I'm seeing in the runtimes (just to get the search a bit more predictable) although that's probably a bit less general and more circuit specific.
I'm not sure if you've discussed any thoughts on this internally?
I also noticed you've started working on an exact search for the routing problem (there's a CVRP cut generator down there) - have you run any performance comparisons between this and the best known approaches?
— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub https://github.com/google/or-tools/issues/912#issuecomment-447360894, or mute the thread https://github.com/notifications/unsubscribe-auth/AKj17fzeKQexpU7KPuOeW9YzR2z6hnveks5u48VHgaJpZM4YBaXJ .
Hope you don't mind me posting this question here but it is related to the circuit constraint (which is central to this discussion).
The SAT layer currently models a directed formulation for the circuit constraint - which is a reasonable way to model the circuit. The other way is to transform the graph to an undirected graph (a somewhat less intuitive representation with n additional dummy edges and nodes) but the graph supports faster operations when doing k-opt style local search.
I'm not sure what the implications are in terms of the propagators that SAT uses - finding literature on this exact topic is a bit tricky. Would/does a transformation such as this impact the SAT search?
It will impact search as we do not have dedicated propagators, not linearization, nor cut generator for the modeling you propose.
This being said, another representation will definitely impact search as you change the nature of the decisions you take. This is unknown territory for us. Laurent Perron | Operations Research | lperron@google.com | (33) 1 42 68 53 00
Le mer. 19 déc. 2018 à 12:40, TheCoolRob notifications@github.com a écrit :
Hope you don't mind me posting this question here but it is related to the circuit constraint (which is central to this discussion).
The SAT layer currently models a directed formulation for the circuit constraint - which is a reasonable way to model the circuit. The other way is to transform the graph to an undirected graph (a somewhat less intuitive representation with n additional dummy edges and nodes) but the graph supports faster operations when doing k-opt style local search.
I'm not sure what the implications are in terms of the propagators that SAT uses - finding literature on this exact topic is a bit tricky. Would/does a transformation such as this impact the SAT search?
— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub https://github.com/google/or-tools/issues/912#issuecomment-448565296, or mute the thread https://github.com/notifications/unsubscribe-auth/AKj17TaBRRnbBZSk4czIe5kNrGD4ic0Aks5u6iW2gaJpZM4YBaXJ .
Thanks for the response Laurent. I realise that a bespoke propagator for this style circuit would need to be added - I suppose my question is more around whether you think there is a benefit or not when handling directed/undirected graphs in SAT?
I was curious as to why directed graphs seem to be the flavour of choice in SAT.
Because the model is with next variables (in CP). This allows dimensions (accumulated values along the path like time, travelled distance, capacity used). Laurent Perron | Operations Research | lperron@google.com | (33) 1 42 68 53 00
Le mer. 19 déc. 2018 à 14:11, TheCoolRob notifications@github.com a écrit :
Thanks for the response Laurent. I realise that a bespoke propagator for this style circuit would need to be added - I suppose my question is more around whether you think there is a benefit or not when handling directed/undirected graphs in SAT?
I was curious as to why directed graphs seem to be the flavour of choice in SAT.
— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub https://github.com/google/or-tools/issues/912#issuecomment-448591928, or mute the thread https://github.com/notifications/unsubscribe-auth/AKj17XNTEvfE9T1jL10AuHXy72tSbh2Zks5u6jrqgaJpZM4YBaXJ .
Heading back to the original topic. I've finished setting up an exact search for the CVRP - the performance is not half bad (about 10 seconds for 50 nodes, 5 vehicles). I've combined the tsp cuts, with about 5 types of capacity cuts.
I picked up a peculiar behaviour though in the results from glop. When the search is backtracking - I would get variable values as the current incumbent that are non-binary. I.e. a boolean variable will come through with a value of 1.5 (as an example) which feels like it's broken the implicit constraint bounds on that column. If I don't attempt to generate any cuts at these LP solutions I'm able to get a valid optimality proof.
I tried clearing the state on the simplex before each call (to check whether that algorithm has an issue) but the behaviour is unchanged. If I assume there is no issue in the configuration of the LP matrix (i.e. they're all configured correctly as [0,1] continuous variables) then does this incumbent represent an infeasible solution? I.e. if I've tightened the bounds on other variables through a cut that causes a particular set of variables to exceed their row constraint on that iteration of the LP? Is the expectation then that the cut generator should verify if the solution is a valid incumbent?
I noticed that in the way that you generate the TSP cuts using the SCC you wouldn't pick this up. Perhaps SAT immediately rejects the solution after the relaxation has been generated - but it feels like there is a missing partial-feasibility/propagation step w.r.t the LP that is missing?
ah! found the culprit. You have a hard-coded limit on the number of simplex iterations.
// TODO(user): Put more at the root, and less afterwards?
parameters.set_max_number_of_iterations(500);
Increasing this solved the problem. I understand the intention here, it's not unreasonable to amortise the simplex convergence over multiple runs, except that you could have infeasible partial solutions which get passed to a cut generator. I think it's just good to know that it's up to the cut generator to then verify the incumbent before working with it.
So I went back and made some changes to the tsp cut generators based on the way glop returns an incumbent. I managed to reduce the search time a fair bit further by adding a higher number of simplex iterations closer to the root node (cuts generated higher in the tree are obviously more useful). I was profiling where the runtime is spent on larger problems (around 150 nodes) and found that out of a 30 second search, only 3 seconds was actually in the search itself and the balance was in the sat::PresolveCpModel (27 seconds). I tried disabling the presolve using
parameters.set_cp_model_use_sat_presolve(false);
but it didn't seem to impact the use of that method. Not sure if you had seen this behaviour before or if it's a know issue with larger matricies with the circuit constraint?
I've been thinking about a small checklist of features which I've seen in other environments that would be really nice to have if you're looking to improve the exact-search on the SAT solver:
Do you already have support for any of these functions or are they already on a roadmap?
Local cuts, not on the roadmap Cut pool manager: being implemented. We are looking at knapsack cuts as a starting point. Column generation: not sure, we may look at auto-lbbd as a way to scale. Not high priority right now. Laurent Perron | Operations Research | lperron@google.com | (33) 1 42 68 53 00
Le mar. 1 janv. 2019 à 13:38, TheCoolRob notifications@github.com a écrit :
So I went back and made some changes to the tsp cut generators based on the way glop returns an incumbent. I managed to reduce the search time a fair bit further by adding a higher number of simplex iterations closer to the root node (cuts generated higher in the tree are obviously more useful). I was profiling where the runtime is spent on larger problems (around 150 nodes) and found that out of a 30 second search, only 3 seconds was actually in the search itself and the balance was in the sat::PresolveCpModel (27 seconds). I tried disabling the presolve using
parameters.set_cp_model_use_sat_presolve(false);
but it didn't seem to impact the use of that method. Not sure if you had seen this behaviour before or if it's a know issue with larger matricies with the circuit constraint?
I've been thinking about a small checklist of features which I've seen in other environments that would be really nice to have if you're looking to improve the exact-search on the SAT solver:
- The current approach is to branch on the most fractional variable when using the LP linearisation - which is a reasonable approach. A more sophisticated approach (on TSPs specifically) is to branch on cliques (so node-based rather than edge based) - is this something that's possible in SAT?
- Could the SAT search support local cuts (all cuts seem to be global at the moment)? This would mean that when moving up the search tree, rows constraints could be dropped from the LP which keeps the simplex running faster.
- Adding a cut pool manager - this is useful for aggregating together cuts, checking for dominant cuts that are found over time (so re-form the LP with far fewer rows), checking for dual representations which are more compact (such as clique inversions) and for duplicate cuts that might be returned from certain algorithms which are already represented in the LP. The best place for this functionality would probably be between the lp-constraint and glop.
- Could SAT support column generation during the search? The use-case for this is a bit more advanced, but on certain problem domains (i.e. dense graphs) you do not need to handle all variables simultaneously and have bounds around groups of variables for the minimum contribution to the objective function (for minimisation) for using any of those variables - you then have a dummy variable in the LP which represents using any one of the variables in this collection (based off the bound) - which means you can dynamically introduce variables into the LP on the fly and change the bound on the dummy variable. This functionality would overlap with the branching strategy for cliques (I suspect one way of implementing this would be to add a new binary variable to the LP - unless SAT has a way already).
Do you already have support for any of these functions or are they already on a roadmap?
— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub https://github.com/google/or-tools/issues/912#issuecomment-450726453, or mute the thread https://github.com/notifications/unsubscribe-auth/AKj17azChpHLdGPsFuOmIIGWQBerUyDGks5u-1bLgaJpZM4YBaXJ .
Thanks for the responses Laurent.
1) No :-) 2) The code is not yet stable, I should publish something this week. 3) https://optimise.amsi.org.au/wp-content/uploads/sites/52/2017/07/stuckey-peter-uom-opt17-talk.pdf
The LBBD looks pretty interesting. I'm familiar with Benders - I quite like the extension of pulling cuts back from the subproblems. Thanks for the reference.
On the clique-style branching: I can hack something in by creating a fixed set of dummy [0,1] variables when the LP initialises and allow it to branch on these variables at specific points in the search (i.e. when a clique branch is available, use it, otherwise default to the most fractional edge) - the variables will simply be unconstrained until the cliques are identified. I know the LP will support adding columns during the search - does SAT allow new variables to be created on the fly during the search? It would be less of a hack if this were possible.
You can add Boolean variables (not in the model proto) during search. I guess IntegerVariable too. Laurent Perron | Operations Research | lperron@google.com | (33) 1 42 68 53 00
Le jeu. 3 janv. 2019 à 08:35, TheCoolRob notifications@github.com a écrit :
The LBBD looks pretty interesting. I'm familiar with Benders - I quite like the extension of pulling cuts back from the subproblems. Thanks for the reference.
On the clique-style branching: I can hack something in by creating a fixed set of dummy [0,1] variables when the LP initialises and allow it to branch on these variables at specific points in the search (i.e. when a clique branch is available, use it, otherwise default to the most fractional edge)
- the variables will simply be unconstrained until the cliques are identified. I know the LP will support adding columns during the search - does SAT allow new variables to be created on the fly during the search? It would be less of a hack if this were possible.
— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub https://github.com/google/or-tools/issues/912#issuecomment-451071977, or mute the thread https://github.com/notifications/unsubscribe-auth/AKj17W35V93YjCFP-UOm9SZ5BNWppJX3ks5u_bK5gaJpZM4YBaXJ .
Hi @lperron. Just popping an update here. I see you've added a ton of new cut-generators in the last year. I wasn't able to re-run my previous benmark because of the changes on the linear_constraint_manager (I was building the lp model directly and it seems that's quite a bit harder with the internal changes now). That said, I was able to find the same problem set I ran previously and pushed it through the cp_model_proto layer and let your cuts do the heavy lifting. It's a serious improvement over what was there before (which I'm sure you know). Your method uses around 40% less time than the approach I glued together above (the red line in the previous chart). Just for comparison:
I'm still keen to have a custom cut-callback in the layer, which I know isn't possible through the proto model - but perhaps if one adds a custom callback to the model before calling the solve that might do the trick? I know presolve is an issue, but perhaps if it it's required that it be disabled if you're going to use an advanced interface like this?
Model model;
model.Add(NewSatParameters(parameters));
model.Add(NewCutCallback(callback1));
const CpSolverResponse response = SolveCpModel(cp_model, &model);
You probably know of a better way from a design perspective? this was just an initial thought.
p.s. I think part of the speed improvement of your current implementation over what I kobbled together is linked to the removing basic variables - I noticed you added that routine for cut removal. It makes a significant difference to the simplex performance.
Hi
Following on from #666 - I found the reason why I wasn't using the
CPModelProto
. I've got some custom cutting planes I'm implementing over a SAT model and using the CutGenerator callbacks.I wasn't able to find an interface on the
CPModelProto
that allows for adding a cut-callback. I would prefer to work at the CPModel level rather than SAT (as per Laurents' suggestion) but perhaps I'm just missing the part of the API that gives the appropriate functionality?It could be that it's not available at that level - passing a pointer to a callback over the proto model seems like a bad idea and something that probably wouldn't be supported. I suppose this issue is simply then for confirmation that this is the case or not.
Thanks