Get multiply feasible solutions from solution pool for mip
AnsweredHi, I am trying to revised tsp_c++.cpp to get multiply feasible solution. I set PoolSearchMode by default, and when I trying to get feasible solution, I just got infeasible solution which has subtour. I am not sure if the callback for subtour elimination was failed or default PoolSearchMode may provide an infeasible solution?
-
Hi,
All solutions in the solution pool should be feasible solutions.
Could you share a minimal reproducible example?
Best, Marika0 -
Hi Marika, Thank you very much for you reply. I just found I could get optimal solution which is definitely feasible when I set PoolSearchMode by default and only get one solution. But, when I want to get more feasible solution with default PoolSearchMode, I found some solution is infeasible, it seems failed to add lazy constraints(by callback). I just don't know why sometimes adding lazy constraints works, sometimes not.(I set GRB_IntParam_LazyConstraints to 1). Since I revised some part of code, I will try to provide an example.
if (len < n) {// cout << "XXXX len = " << len << endl;// for (int i = 0; i < len; i++)// cout << tour[i] << " ";// cout << endl;// for (int i = 0; i < n; i++){// for (int j = i+1; j < n; j++) {// if (x[i][j] > 0.5)// cout << "x:" << i << " " << j << endl;// if (x_bar[i][j] > 0.5)// cout << "x_bar:" << i << " " << j << endl;// }// }// Add subtour elimination constraintGRBLinExpr expr = 0;//cout << "starting building constraints.." << endl;for(int i =0; i < len-1; i++){for(int j = i+1; j < len; j++){expr += vars[tour[i]][tour[j]];expr += vars_bar[tour[i]][tour[j]];}}//cout << "finishing building constraints.." << endl;addLazy(expr <= len-1);//cout << "finishing add lazy constraints.." << endl;}0 -
This is not a reproducible example.
If I understand you correctly you started with the Gurobi example tsp_c++.cpp. But why do you change something in the code if you want to get more feasible solutions with default PoolSearchMode?
If I run the example tsp_c++.cpp with value 12, Gurobi finds 2 solutions:
Solution count 2: 3.41529 6.06623
Both solutions are feasible. If you want Gurobi to find more solutions, you need to change the parameter PoolSearchMode to 1 or 2 (and maybe also adapt PoolSolutions).
I do not see where you get infeasible solutions.0 -
Yes. I could also get the results you shown here. But if I want to print these two solutions with GRB_IntParam_SolutionNumber, sometimes the sub-optimal solution is infeasible, which means it has subtour, like "0-1-3-2" and "4-6-7-5" if we have 8 nodes. Do you mean I have to set PoolSearchMode to 1 or 2 to get feasible sub-optimal solution?
0 -
Since I need to revised tsp_c++ for my own project, I just show you the part of code and infeasible result I got. For example, I have 22 node, but I got a tour with only 5 nodes.
int SolCount = model.get(GRB_IntAttr_SolCount);
if(SolCount >0){
double poolobj =0.0;
for(int i =0; i < SolCount; i++){
model.set(GRB_IntParam_SolutionNumber, i);
poolobj = model.get(GRB_DoubleAttr_PoolObjVal);
double**sol =newdouble*[n];
double**sol_bar =newdouble*[n];
for(int i =0; i < n; i++)
{
sol[i]=model.get(GRB_DoubleAttr_Xn,vars[i], n);
sol_bar[i]=model.get(GRB_DoubleAttr_Xn,vars_bar[i], n);
}
int* tour =newint[n];
int len;
findsubtour(n, sol, sol_bar,&len, tour);
if(len < n){
cout << "len = " << len << endl;
cout << "tour : " << " ";
for(int i =0; i < len; i++){
cout << tour[i] << " ";
}
cout << endl;The result is as following.
len = 5
tour : 3 4 11 16 19
sol:0 14
sol:0 17
sol:1 2
sol:1 6
sol:2 5
sol:3 4
sol_bar:3 19
sol:4 11
sol:5 7
sol:6 10
sol:7 9
sol:8 10
sol_bar:8 13
sol_bar:9 17
sol_bar:11 16
sol_bar:12 13
sol:12 15
sol:14 21
sol:15 18
sol:16 19
sol:18 20
sol:20 210 -
This is still not enough information to help you. Your problem is probably not defined properly.
Could you explain a bit more about what you exactly want to solve?
Why do you have variables vars and vars_bar? What does vars_bar represent?
What does the callback to add the sub-tour elimination constraints look like?
What does the output "sol: 0 14" or "sol_bar: 3 19" mean?0 -
Thank you very much for your reply! Since I revised tsp_c++ to fit my model, vars and vars_bar are both decision variables. For example, as for 0 14, vars[0][14] + vars_bar[0][14] = 1, which means if we travel from 0 to 14, we could make vars[0][14] = 1 or vars_bar[0][14] = 1. As for the solution part, sol: 0 14 means the optimal solution told us that we travel from 0 to 14 by vars[0][14] = 1 and we travel from 3 to 19 by vars_bar[3][19] = 1. From above result, when I get this suboptimal solution, this suboptimal solution has a subtour 3 4 11 16 19, which is infeasible because subtour solution should be eliminated by lazy constraints. So I think the callback function didn't work for this suboptimal, but I have no idea how to fix it. Here is the callback part. Thank you.
void callback() {
try{
if(where == GRB_CB_MIPSOL){
// Found an integer feasible solution - does it visit every node?
double**x =newdouble*[n];
double**x_bar =newdouble*[n];
int*tour =newint[n];
int len;
for(int i =0; i < n; i++){
x[i]=getSolution(vars[i], n);
x_bar[i]=getSolution(vars_bar[i],n);
}
findsubtour(n, x, x_bar,&len, tour);
if(len < n){
// Add subtour elimination constraint
GRBLinExpr expr = 0;
for(int i =0; i < len-1; i++){
for(int j = i+1; j < len; j++){
expr += vars[tour[i]][tour[j]];
expr += vars_bar[tour[i]][tour[j]];
}
}
addLazy(expr <= len-1);
}
for(int i =0; i < n; i++)
{
delete[]x[i];
delete[]x_bar[i];
}
delete[] x;
delete[] x_bar;
delete[] tour;
}
}catch(GRBException e){
cout << "Error number: " << e.getErrorCode() << endl;
cout << e.getMessage() << endl;
}catch(...){
cout << "Error during callback" << endl;
}
}0 -
I cannot reproduce the issue.
You probably need to share the complete reproducible example.
0 -
Hi, Here is my code. The main function of this project is to find several routes that start from 0 and visited several nodes and finally go back to 0. For example, if I have 0 1 2 3 4 5 6 7 8, and I would like to get the results that looks like 0-1-2-4-5-0 and 0-3-8-7-6-0. And I do not want to get subtour like 3-8-7-3. When I retrieved the suboptimal solution. I did found certain suboptimal solution including subtour like 3-8-7-3. Please see the code and set node number 33.
#include "gurobi_c++.h"
And this is one of suboptimal solution I got from solution pool that has sub tour.
#include <cassert>
#include <cstdlib>
#include <cmath>
#include <sstream>
#include <vector>
#include <algorithm>
using namespace std;
string itos(int i) {stringstream s; s << i; return s.str(); }
double distance(double* x, double* y, int i, int j);
void findsubtour(int n, double** sol, vector<vector<int>> &sub_tour, int* tourlenP);
// Subtour elimination callback. Whenever a feasible solution is found,
// find the smallest subtour, and add a subtour elimination constraint
// if the tour doesn't visit every node.
class subtourelim: public GRBCallback
{
public:
GRBVar** vars;
intn;
subtourelim(GRBVar**xvars,intxn){
vars=xvars;
n=xn;
}
protected:
voidcallback(){
try{
if(where == GRB_CB_MIPSOL){
// Found an integer feasible solution - does it visit every node?
double**x=newdouble*[n];
for(inti=0;i<n;i++)
x[i]=getSolution(vars[i],n);
inttour_len;
vector<vector<int>>sub_tour;
findsubtour(n,x,sub_tour,&tour_len);
for(size_ti=0;i<sub_tour.size();i++){
if(find(sub_tour[i].begin(),sub_tour[i].end(),0)==sub_tour[i].end()){
intlen=sub_tour[i].size();
GRBLinExpr expr = 0;
for(size_tk=0;k<sub_tour[i].size();k++)
for(size_tj=k+1;j<sub_tour[i].size();j++)
expr+=vars[sub_tour[i][k]][sub_tour[i][j]];
addLazy(expr<=len-1);
cout << "Adding lazy constraint" << endl;
}
}
for(inti=0;i<n;i++)
delete[]x[i];
delete[]x;
}
}catch(GRBException e){
cout << "Error number: " << e.getErrorCode() << endl;
cout << e.getMessage() << endl;
}catch(...){
cout << "Error during callback" << endl;
}
}
};
// Given an integer-feasible solution 'sol', find the smallest
// sub-tour. Result is returned in 'tour', and length is
// returned in 'tourlenP'.
void
findsubtour(int n,
double**sol,
vector<vector<int>>&sub_tour,
int*tourlenP)
{
bool*seen=newbool[n];
for(inti=0;i<n;i++)
seen[i]=false;
intvisited=0;
boolflag;
vector<int>temp_tour;
do{
flag=false;
intnode=0;
seen[node]=true;
for(intlen=0;len<n;len++){
inti;
for(i=0;i<n;i++){
if(sol[node][i]>0.5&&!seen[i]){
if(!node){
temp_tour.emplace_back(node);
}
node=i;
temp_tour.emplace_back(node);
seen[node]=true;
++visited;
break;
}
}
if(i==n){
sub_tour.emplace_back(temp_tour);
temp_tour.clear();
break;
}
}
for(intj=1;j<n;j++){
if(sol[0][j]>0.5&&!seen[j]){
flag=true;
break;
}
}
}while(flag);
cout << "visited = " << visited << endl;
*tourlenP=visited;
while(visited<n-1){
for(inti=1;i<n;i++){
while(!seen[i]){
intstart=i;
for(intlen=0;len<n;len++){
temp_tour.emplace_back(start);
++visited;
seen[start]=true;
intj;
for(j=1;j<n;j++){
if(sol[start][j]>0.5&&!seen[j]){
start=j;
break;
}
}
if(j==n){
sub_tour.emplace_back(temp_tour);
temp_tour.clear();
break;
}
}
}
}
}
// cout << "findsubtour : sub_tour_list" << endl;
// for (size_t i = 0; i < sub_tour.size(); i++){
// for (size_t j = 0; j < sub_tour[i].size(); j++){
// cout << sub_tour[i][j] << " ";
// }
// cout << endl;
// }
delete[]seen;
}
// Euclidean distance between points 'i' and 'j'.
double
distance(double* x,
double*y,
inti,
intj)
{
doubledx=x[i]-x[j];
doubledy=y[i]-y[j];
returnsqrt(dx*dx+dy*dy);
}
int
main(int argc,
char*argv[])
{
if(argc<2){
cout << "Usage: tsp_c++ size" << endl;
return1;
}
intn=atoi(argv[1]);
double*x=newdouble[n];
double*y=newdouble[n];
inti;
for(i=0;i<n;i++){
x[i]=((double)rand())/RAND_MAX;
y[i]=((double)rand())/RAND_MAX;
}
GRBEnv *env = NULL;
GRBVar **vars = NULL;
vars = new GRBVar*[n];
for(i=0;i<n;i++)
vars[i]=newGRBVar[n];
try{
intj;
env = new GRBEnv();
GRBModel model = GRBModel(*env);
// Must set LazyConstraints parameter when using lazy constraints
model.set(GRB_IntParam_LazyConstraints,1);
//model.set(GRB_IntParam_PoolSearchMode,1);
// Create binary decision variables
for(i=0;i<n;i++){
for(j=0;j<=i;j++){
vars[i][j]=model.addVar(0.0,1.0,distance(x,y,i,j),
GRB_BINARY, "x_"+itos(i)+"_"+itos(j));
vars[j][i]=vars[i][j];
}
}
// Degree-2 constraints
GRBLinExpr expr1 = 0;
for(j=0;j<n;j++){
expr1+=vars[0][j];
}
model.addConstr(expr1>=4,"vehicle_");
for(i=1;i<n;i++){
GRBLinExpr expr = 0;
for(j=0;j<n;j++)
expr+=vars[i][j];
model.addConstr(expr==2,"deg2_"+itos(i));
}
// Forbid edge from node back to itself
for(i=0;i<n;i++)
vars[i][i].set(GRB_DoubleAttr_UB,0);
// Set callback function
subtourelimcb=subtourelim(vars,n);
model.setCallback(&cb);
// Optimize model
model.optimize();
// model.write("tsp_model.lp");
cout << "########Print the results#########" << endl;
// Extract solution
intSolCount=model.get(GRB_IntAttr_SolCount);
if(SolCount>0){
for(inti=0;i<SolCount;i++){
model.set(GRB_IntParam_SolutionNumber,i);
doublepoolobj=model.get(GRB_DoubleAttr_PoolObjVal);
cout << "Obj value " << i << ": " << poolobj << " ";
double**sol=newdouble*[n];
for(inti=0;i<n;i++)
sol[i]=model.get(GRB_DoubleAttr_Xn,vars[i],n);
// for (int i = 0; i < n; i++){
// for (int j = i+1; j < n; j++) {
// if (sol[i][j] > 0.5)
// cout << i << " " << j << " ";
// }
// }
// cout << endl;
intvisited_len;
vector<vector<int>>sol_tour;
findsubtour(n,sol,sol_tour,&visited_len);
cout << "visited_len = " << visited_len << endl;
//assert(len == n);
cout << "sub_tour_list" << endl;
for(size_ti=0;i<sol_tour.size();i++){
for(size_tj=0;j<sol_tour[i].size();j++){
cout << sol_tour[i][j] << " ";
}
cout << endl;
}
for(inti=0;i<n;i++)
delete[]sol[i];
delete[]sol;
}
}
}catch(GRBException e){
cout << "Error number: " << e.getErrorCode() << endl;
cout << e.getMessage() << endl;
}catch(...){
cout << "Error during optimization" << endl;
}
for(i=0;i<n;i++)
delete[]vars[i];
delete[] vars;
delete[]x;
delete[]y;
delete env;
return0;
}sub_tour_list
0 17 8 5 6 22
0 20 15 18 21 3 19 30 24 11 26 9 4 32 12 10 28 13 29 25
0 23 1 7 31
2 14 27 160 -
Which Gurobi Version do you use?
Could you please also add the log file?
0 -
The Gurobi version I use is 952. When you run this code, you are able to get the log file.
0 -
There were some fixes for Version 10.0.0, see Recent Bug Fixes by Version - Gurobi Optimization:
- Fixed issue with not calling the callback for solutions found in MIP cleanup at the end of the solve
- Fixed a bug that the MIPSOL callback is not called for one of the solutions found
I think the latter corresponds to what you observed.
I recommend updating to the newest version.
0 -
Thank you very much! Problem solved!
0
Please sign in to leave a comment.
Comments
13 comments