qbpp_grb.hpp

Go to the documentation of this file.

 
#ifndef GUROBI_QBPP_HPP
#define GUROBI_QBPP_HPP
 
#include <cmath>
 
#include "abs2.hpp"
#include "gurobi_c++.h"
#include "qbpp.hpp"
 
#define GRB_SAFE_CALL(func)                                               \
  try {                                                                   \
    func;                                                                 \
  } catch (GRBException e) {                                              \
    std::cerr << e.getErrorCode() << ": " << e.getMessage() << std::endl; \
    exit(1);                                                              \
  }
 
namespace qbpp_grb {
 
//===========================
// Class forward declarations
//===========================
class QuadModel;
class Sol;
class Callback;
 
//===================
// Class declarations
//===================
 
class QuadModel : public qbpp::QuadModel {
 protected:
  GRBEnv grb_env;
  std::shared_ptr<GRBModel> grb_model_ptr;
  GRBVar *grb_x;
 
 public:
 
  QuadModel(const qbpp::QuadModel &quad_model, bool verbose = false);
 
  QuadModel(const QuadModel &grb_model) = default;
 
  void set(const std::string &key, const std::string &val) {
    grb_model_ptr->getEnv().set(key, val);
  };
 
  void set_time_limit(uint32_t time_limit) {
    set("TimeLimit", std::to_string(time_limit));
  }
 
  void set(Callback &cb);
 
 
  GRBVar get_grb_var(int index) const { return grb_x[index]; };
 
  GRBModel &get_grb_model() { return *grb_model_ptr; }
 
  const GRBModel &get_grb_model() const { return *grb_model_ptr; }
 
  Sol optimize();
 
  void write(std::string filename) {
    GRB_SAFE_CALL(grb_model_ptr->write(filename));
  }
};
 
class Sol : public qbpp::Sol {
  qbpp::energy_t bound;
 
 public:
 
  Sol(const qbpp::QuadModel &quad_model) : qbpp::Sol(quad_model) {}
 
 
  qbpp::energy_t get_bound() const { return bound; }
 
  qbpp::energy_t set_bound(qbpp::energy_t eval) {
    bound = eval;
    return bound;
  }
};
 
class Callback : public GRBCallback {
 protected:
  const QuadModel quad_model;
 
  const GRBModel &grb_model;
 
  std::optional<qbpp::energy_t> target_energy = std::nullopt;
 
  //  std::chrono::high_resolution_clock::time_point start_time =
  //      std::chrono::high_resolution_clock::now();
 
  mutable std::mutex mtx;
 
 public:
  Callback(const QuadModel &quad_model)
      : quad_model(quad_model), grb_model(quad_model.get_grb_model()) {
    qbpp::get_time();
  };
 
  virtual ~Callback() = default;
 
  Sol get_sol();
 
  virtual void callback() override {
    if (where == GRB_CB_MIPSOL) {
      qbpp::energy_t energy = get_sol().get_energy();
      std::cout << "TTS = " << std::fixed << std::setprecision(3)
                << std::setfill('0') << qbpp::get_time()
                << "s Energy = " << energy << std::endl;
      abort_if_target_energy(energy);
    }
  }
 
  void set_target_energy(qbpp::energy_t target_energy) {
    this->target_energy = target_energy;
  }
 
  void abort_if_target_energy(qbpp::energy_t energy) {
    if (target_energy.has_value() && energy <= target_energy) {
      abort();
    }
  }
 
  double getDoubleInfoPublic(int what) { return getDoubleInfo(what); }
 
  double getSolutionPublic(GRBVar v) { return getSolution(v); }
};
 
//=============================
// Class QuadModel member functions
//=============================
 
inline QuadModel::QuadModel(const qbpp::QuadModel &quad_model, bool verbose)
    : qbpp::QuadModel(quad_model), grb_env(true) {
  // Suppress Gurobi outputs to the screen.
  grb_env.set(GRB_IntParam_OutputFlag, 0);
  if (verbose) {
    grb_env.set("OutputFlag", "1");
  }
  GRB_SAFE_CALL(grb_env.start());
 
  grb_model_ptr = std::make_unique<GRBModel>(grb_env);
 
  //  Objective is minimization.
  GRB_SAFE_CALL(grb_model_ptr->set(GRB_IntAttr_ModelSense, GRB_MINIMIZE));
 
  // Add binary variables to the Gurobi model.
  grb_x = grb_model_ptr->addVars(var_count(), GRB_BINARY);
  // Compute objective function of QUBO problem.
  GRBQuadExpr obj;
  obj += quad_model.get_constant();
  for (qbpp::vindex_t i = 0; i < var_count(); ++i) {
    if (quad_model.get_linear(i) != 0)
      obj += quad_model.get_linear(i) * grb_x[i];
  }
  for (qbpp::vindex_t i = 0; i < var_count(); ++i) {
    for (qbpp::vindex_t j = 0; j < quad_model.get_degree(i); ++j) {
      auto [k, coeff] = quad_model.get_quadratic(i, j);
      if (i < k) obj += coeff * grb_x[i] * grb_x[k];
    }
  }
 
  // Set the objective function to the Gurobi model.
  GRB_SAFE_CALL(grb_model_ptr->setObjective(obj));
}
 
// void QuadModel::set(Callback *cb) { grb_model_ptr->setCallback(cb); }
 
inline void QuadModel::set(Callback &cb) { grb_model_ptr->setCallback(&cb); }
 
inline Sol QuadModel::optimize() {
  Sol sol(*this);
  GRB_SAFE_CALL(grb_model_ptr->optimize());
  sol.set_energy(std::round(grb_model_ptr->get(GRB_DoubleAttr_ObjVal)));
  sol.set_bound(std::round(grb_model_ptr->get(GRB_DoubleAttr_ObjBound)));
  for (qbpp::vindex_t i = 0; i < var_count(); ++i)
    sol.set(i, int(grb_x[i].get(GRB_DoubleAttr_X)));
  return sol;
}
 
//================================
// Class Callback member functions
//================================
 
inline Sol Callback::get_sol() {
  Sol sol(quad_model);
  std::lock_guard<std::mutex> lock(mtx);
  // Intentionally commented out because the solution obtained by getSolution
  // may provide corrupted values when Gurobi finds an optimal solution.
  // sol.set_energy(std::round(getDoubleInfoPublic(GRB_CB_MIPSOL_OBJ)));
  sol.set_bound(std::round(getDoubleInfoPublic(GRB_CB_MIPSOL_OBJBND)));
  // Get the solution from Gurobi Optimizer.
  for (qbpp::vindex_t i = 0; i < quad_model.var_count(); ++i) {
    sol.set(i, int(getSolutionPublic(quad_model.get_grb_var(i))));
  }
  // eval is used to compute the energy of the solution instead of MIPSOL_OBJ.
  sol.set_energy(eval(quad_model, sol));
  return sol;
}
 
}  // namespace qbpp_grb
#endif  // __GUROBI_QBPP_HPP__