CGAL / cgal

The public CGAL repository, see the README below
https://github.com/CGAL/cgal#readme
Other
5.03k stars 1.39k forks source link

disable presolver when trying to solve mixed_integer_program_traits with GLPK #5494

Closed maximouth closed 3 years ago

maximouth commented 3 years ago

Please use the following template to help us solving your issue.

Issue Details

Hello,

I'm trying to solve a MILP problem with CGAL using the GLPK_mixed_integer_program_traits class. After defining my problem, all the variables, the different constraints and the linear objective function I launch the solver. But as result I got the following sentence : " solving problem failed error: Unable to start the search, because LP relaxation of the MIP problem instance hasno primal feasible solution. (This code may appear only if the presolver is enabled.)"

Can you tell me how to disable the presolver for the resolution of this problem ?

Thank you in advance, Maxime

Source Code

#include <iostream>
#ifdef CGAL_USE_SCIP
#include <CGAL/SCIP_mixed_integer_program_traits.h>
typedef CGAL::SCIP_mixed_integer_program_traits<double>                        MIP_Solver;
#elif defined(CGAL_USE_GLPK)
#include <CGAL/GLPK_mixed_integer_program_traits.h>
typedef CGAL::GLPK_mixed_integer_program_traits<double>                        MIP_Solver;
#endif
#if defined(CGAL_USE_GLPK) || defined(CGAL_USE_SCIP)
typedef typename MIP_Solver::Variable                        Variable;
typedef typename MIP_Solver::Linear_objective        Linear_objective;
typedef typename MIP_Solver::Linear_constraint        Linear_constraint;

int main(int argc, char* argv[])
{
  MIP_Solver solver;

  if (argc != 2) {
    std::cout << "Wrong usage, 1 arg required :\n ./dobss \"Number of nodes\"\n"; 
    return -1;
  }

  // The number of node to defend inside the airport  
  int node_number = atoi(argv[1]);

  double k = node_number;

  double M = 123456;

  //The reward matrix for the defender
  double Rij [node_number][node_number];
  //The reward matrix for the attacker
  double Cij [node_number][node_number];    

  // Zi_j variable
  Variable *zij[node_number][node_number];

  Variable *Qj[node_number];
  Variable *A;
  int name;

  for (int i = 0 ; i < node_number ; i++) {
    for (int j = 0 ; j < node_number ; j++) {
      zij[i][j] = solver.create_variable(Variable::INTEGER, 0., node_number, "Z" + std::to_string(i) + std::to_string(j));
    }
  }

  for (int i = 0 ; i < node_number ; i++) {
    Qj[i] = solver.create_variable(Variable::BINARY, 0, 1, "Q" + std::to_string(i));
  }

  A = solver.create_variable(Variable::CONTINUOUS, -Linear_constraint::infinity(), Linear_constraint::infinity(), "A");

  for (int i = 0 ; i < node_number ; i++) {
    for (int j = 0 ; j < node_number ; j++) {
      // When the attacker and the defender target the same node
      // The defender wins 1
      // The attacker looses 1
      if (i == j) {
    Rij[i][j] = 1.0;
    Cij[i][j] = -1.0;
      }
      // When the attacker and the defender don't target the same node
      // The defender looses 1
      // The attacker wins 1
      else {
    Rij[i][j] = -1.0;
    Cij[i][j] = 1.0;    
      }
    }
  }

  // Objective.
  // Be careful this is "MAXIMIZE"
  Linear_objective * obj = solver.create_objective(Linear_objective::MAXIMIZE);

for (int i = 0 ; i < node_number ; i++) {
 obj->add_coefficient(Qj[i], 0.0);
  for (int j = 0 ; j < node_number ; j++) {
    obj->add_coefficient(zij[i][j], (1.0/node_number) * Rij[i][j]);
  }
 }

 obj->add_coefficient(A, 0.0);

 // Constraint c1: \forall i \forall j : \sum Zij = k
 Linear_constraint* c1 = solver.create_constraint(node_number, node_number, "c1");

 for (int i = 0 ; i < node_number ; i++) {
   c1->add_coefficient(Qj[i], 0.0);
   for (int j = 0 ; j < node_number ; j++) {
     c1->add_coefficient(zij[i][j], 1.0);
   }
 }

 c1->add_coefficient(A, 0.0);

 // Constraint c2: \forall j \sum Zij \leq k
 Linear_constraint* c2;

 for (int k = 0 ; k < node_number ; k++) {
   c2 = solver.create_constraint(-Linear_constraint::infinity(), node_number, "c" + std::to_string(k+2));

   for (int i = 0 ; i < node_number ; i++) {
     c2->add_coefficient(Qj[i], 0.0);
     for (int j = 0 ; j < node_number ; j++) {
       if (i == k) {
     c2->add_coefficient(zij[i][j], 1.0);
       }
       else {
     c2->add_coefficient(zij[i][j], 0.0);
       }
     }
   }
   c2->add_coefficient(A, 0.0);
 }

  // Constraint c3: k*q_j \leq \forall i \sum Zij 
 Linear_constraint* c3;

 for (int k = 0 ; k < node_number ; k++) {
   c3 = solver.create_constraint(0, Linear_constraint::infinity(), "c" + std::to_string(k+2+(node_number * node_number)));

   for (int i = 0 ; i < node_number ; i++) {
     c3->add_coefficient(Qj[i], -node_number);
     for (int j = 0 ; j < node_number ; j++) {
       if (j== k) {
     c3->add_coefficient(zij[i][j], 1.0);
       }
       else {
     c3->add_coefficient(zij[i][j], 0.0);
       }
     }
   }
   c3->add_coefficient(A, 0.0);
 }

 // Constraint c4: \forall i \sum Zij \leq k 
 Linear_constraint* c4;

 for (int k = 0 ; k < node_number ; k++) {
   c4 = solver.create_constraint(-Linear_constraint::infinity(), node_number , "c" + std::to_string(k+2+(node_number * node_number)));

   for (int i = 0 ; i < node_number ; i++) {
     c4->add_coefficient(Qj[i], 0.0);
     for (int j = 0 ; j < node_number ; j++) {
       if (j== k) {
     c4->add_coefficient(zij[i][j], 1.0);
       }
       else {
     c4->add_coefficient(zij[i][j], 0.0);
       }
     }
   }
   c4->add_coefficient(A, 0.0);
 }

 // Constraint c5: \forall j \sum Qj = 1 
  Linear_constraint* c5;
 c5 = solver.create_constraint(1., 1. , "c" + std::to_string(k+2+(node_number * node_number)+(node_number * node_number)));

  for (int i = 0 ; i < node_number ; i++) {

    c5->add_coefficient(Qj[i], 1.0);
    for (int j = 0 ; j < node_number ; j++) {
      c5->add_coefficient(zij[i][j], 0.0);
    }
  }
  c5->add_coefficient(A, 0.0);

 // Constraint c6: a - \forall i \sum 1/k Cij * \forall h \sum Zij \geq 0
 Linear_constraint* c6;

 for (int k = 0 ; k < node_number ; k++) {
   c6 = solver.create_constraint(0, Linear_constraint::infinity(), "c" + std::to_string(k+2+(node_number * node_number)));

   for (int i = 0 ; i < node_number ; i++) {
     c6->add_coefficient(Qj[i], 0.0);
     for (int j = 0 ; j < node_number ; j++) {
       if (j== k) {
     c6->add_coefficient(zij[i][j], (-1.0/node_number) * Cij[i][j]);
       }
       else {
     c6->add_coefficient(zij[i][j], 0.0);
       }
     }
   }
   c6->add_coefficient(A, 1.0);
 }

 // Constraint c7: a a - \forall i \sum 1/k Cij * \forall h \sum Zij \leq (A- qj) M
 Linear_constraint* c7;

 for (int k = 0 ; k < node_number ; k++) {
   c7 = solver.create_constraint(-Linear_constraint::infinity(), M , "c" + std::to_string(k+2+(node_number * node_number)));

   for (int i = 0 ; i < node_number ; i++) {
     c7->add_coefficient(Qj[i], M);
     for (int j = 0 ; j < node_number ; j++) {
       if (j== k) {
     c7->add_coefficient(zij[i][j], (-1.0/node_number) * Cij[i][j]);
       }
       else {
     c7->add_coefficient(zij[i][j], 0.0);
       }
     }
   }
   c7->add_coefficient(A, 1.0);
 }

 //Print problem info :

 std::cout << "problem info : \n" ;
 std::cout << "Mixed integer problem : " ;
 std::cout << solver.is_mixed_integer_program() << std::endl;

 std::cout << "integer program :";
 std::cout << solver.is_integer_program() << std::endl;

 std::cout << "binary problem : " ;
 std::cout << solver.is_binary_program() << std::endl << std::endl;

 std::cout << "number continuous variable : " ;
 std::cout << solver.number_of_continuous_variables() << std::endl;

 std::cout << "number integer variable : " ;
 std::cout << solver.number_of_integer_variables() << std::endl;

 std::cout << "number binary variable : " ;
 std::cout << solver.number_of_binary_variables() << std::endl << std::endl;

 std::cout << "objective sense : Min=0 ; Max=1 ; UDF : 2 : " ;
 Linear_objective* ol = solver.objective();  
 std::cout << ol->sense() << std::endl << std::endl;

 std::cout << "objective function : " ;
  for (int i = 0 ; i < solver.number_of_continuous_variables() + solver.number_of_integer_variables() + solver.number_of_binary_variables() ; i++ ) {
   std::cout << solver.variables()[i]->name() << "*" <<ol->get_coefficient(solver.variables()[i]) << " + " ;
 }
 std::cout << std::endl << std::endl;

 std::cout << "constraints : \n" ;
 const std::vector<Linear_constraint*>& cs = solver.constraints();  

 int j = 0;
 while (cs[j] != NULL) {
   std::cout << cs[j]->name() << " lb : " << cs[j]->lower_bound() << ", ub : " << cs[j]->upper_bound() << ". coeff : \n";

   for (int i = 0 ; i < solver.number_of_continuous_variables() + solver.number_of_integer_variables() + solver.number_of_binary_variables() ; i++ ) {
     std::cout << cs[j]->get_coefficient(solver.variables()[i]) << "*"<< solver.variables()[i]->name() << " + ";
   }
   std::cout << std::endl;
   j++;
 }

 std::cout << std::endl;

 // Solve
  if (solver.solve()) {
    std::cout << "result: " << std::endl;
    const std::vector<double>& results = solver.solution();

    for (std::size_t i = 0; i < results.size(); ++i) {
      //std::cout << "\tx" << i + 1 << ": " << results[i] << std::endl;
      std::cout << solver.variables()[i]->name() << " : " << results[i] << std::endl;
    }
    return EXIT_SUCCESS;
  }
  else {
    std::cerr << "solving problem failed" << std::endl;
    std::cout << "error: " << solver.error_message() << std::endl;
    return EXIT_FAILURE;
  }
}
#else
int main(int , char**)
{
  std::cerr << "This test requires either GLPK or SCIP.\n";
  return EXIT_SUCCESS;
}
#endif  // defined(CGAL_USE_GLPK) || defined(CGAL_USE_SCIP)

Environment

sgiraudot commented 3 years ago

Is there a particular reason why you want to use CGAL for this problem? The CGAL GLPK class you are using really just is a wrapper that we only use for one algorithm (Polygonal reconstruction). I'm not sure there's any added value to use this wrapper alone.

Why not directly use GLPK (and having thus access to the lower level functions that you need)?

maximouth commented 3 years ago

Because I didn't knew that this part of CGAL was only a wrapper of GLPK... For this case, using CGAL really adds nothing ? If it is the case, then I will move to GLPK only.

Thanks for the response!