new protobuf to store linear models and solutions, much smaller and faster

makefiles/Makefile.cpp.mk

@@ -474,7 +474,9 @@ LINEAR_SOLVER_LIB_OBJS = \
 	$(OBJ_DIR)/gurobi_interface.$O \
 	$(OBJ_DIR)/linear_solver.$O \
 	$(OBJ_DIR)/linear_solver.pb.$O \
+	$(OBJ_DIR)/linear_solver2.pb.$O \
 	$(OBJ_DIR)/model_exporter.$O \
+        $(OBJ_DIR)/proto_tools.$O \
 	$(OBJ_DIR)/scip_interface.$O \
 	$(OBJ_DIR)/sulum_interface.$O
 
@@ -491,7 +493,7 @@ $(OBJ_DIR)/glpk_interface.$O:$(SRC_DIR)/linear_solver/glpk_interface.cc
 $(OBJ_DIR)/gurobi_interface.$O:$(SRC_DIR)/linear_solver/gurobi_interface.cc
 	$(CCC) $(CFLAGS) -c $(SRC_DIR)/linear_solver/gurobi_interface.cc $(OBJ_OUT)$(OBJ_DIR)$Sgurobi_interface.$O
 
-$(OBJ_DIR)/linear_solver.$O:$(SRC_DIR)/linear_solver/linear_solver.cc $(GEN_DIR)/linear_solver/linear_solver.pb.h
+$(OBJ_DIR)/linear_solver.$O:$(SRC_DIR)/linear_solver/linear_solver.cc $(GEN_DIR)/linear_solver/linear_solver.pb.h $(GEN_DIR)/linear_solver/linear_solver2.pb.h
 	$(CCC) $(CFLAGS) -c $(SRC_DIR)/linear_solver/linear_solver.cc $(OBJ_OUT)$(OBJ_DIR)$Slinear_solver.$O
 
 $(OBJ_DIR)/linear_solver.pb.$O:$(GEN_DIR)/linear_solver/linear_solver.pb.cc
@@ -502,9 +504,20 @@ $(GEN_DIR)/linear_solver/linear_solver.pb.cc:$(SRC_DIR)/linear_solver/linear_sol
 
 $(GEN_DIR)/linear_solver/linear_solver.pb.h:$(GEN_DIR)/linear_solver/linear_solver.pb.cc
 
+$(OBJ_DIR)/linear_solver2.pb.$O:$(GEN_DIR)/linear_solver/linear_solver2.pb.cc
+	$(CCC) $(CFLAGS) -c $(GEN_DIR)/linear_solver/linear_solver2.pb.cc $(OBJ_OUT)$(OBJ_DIR)$Slinear_solver2.pb.$O
+
+$(GEN_DIR)/linear_solver/linear_solver2.pb.cc:$(SRC_DIR)/linear_solver/linear_solver2.proto
+	$(PROTOBUF_DIR)/bin/protoc --proto_path=$(INC_DIR) --cpp_out=$(GEN_DIR) $(SRC_DIR)/linear_solver/linear_solver2.proto
+
+$(GEN_DIR)/linear_solver/linear_solver2.pb.h:$(GEN_DIR)/linear_solver/linear_solver2.pb.cc
+
 $(OBJ_DIR)/model_exporter.$O:$(SRC_DIR)/linear_solver/model_exporter.cc $(GEN_DIR)/linear_solver/linear_solver.pb.h
 	$(CCC) $(CFLAGS) -c $(SRC_DIR)/linear_solver/model_exporter.cc $(OBJ_OUT)$(OBJ_DIR)$Smodel_exporter.$O
 
+$(OBJ_DIR)/proto_tools.$O:$(SRC_DIR)/linear_solver/proto_tools.cc $(GEN_DIR)/linear_solver/linear_solver.pb.h $(GEN_DIR)/linear_solver/linear_solver2.pb.h
+	$(CCC) $(CFLAGS) -c $(SRC_DIR)/linear_solver/proto_tools.cc $(OBJ_OUT)$(OBJ_DIR)$Sproto_tools.$O
+
 $(OBJ_DIR)/scip_interface.$O:$(SRC_DIR)/linear_solver/scip_interface.cc
 	$(CCC) $(CFLAGS) -c $(SRC_DIR)/linear_solver/scip_interface.cc $(OBJ_OUT)$(OBJ_DIR)$Sscip_interface.$O
 

src/linear_solver/glpk_interface.cc

@@ -35,18 +35,17 @@
 #if defined(USE_GLPK)
 
 extern "C" {
-  #include "glpk.h"
+#include "glpk.h"
 }
 
 DECLARE_double(solver_timeout_in_seconds);
 DECLARE_string(solver_write_model);
 
 namespace operations_research {
-
 // Class to store information gathered in the callback
 class GLPKInformation {
  public:
-  explicit GLPKInformation(bool maximize): num_all_nodes_(0) {
+  explicit GLPKInformation(bool maximize) : num_all_nodes_(0) {
     ResetBestObjectiveBound(maximize);
   }
   void Reset(bool maximize) {
@@ -64,8 +63,6 @@ class GLPKInformation {
   double best_objective_bound_;
 };
 
-
-
 // Function to be called in the GLPK callback
 void GLPKGatherInformationCallback(glp_tree* tree, void* info) {
   CHECK_NOTNULL(tree);
@@ -122,8 +119,7 @@ class GLPKInterface : public MPSolverInterface {
   // Change a coefficient in a constraint.
   virtual void SetCoefficient(MPConstraint* const constraint,
                               const MPVariable* const variable,
-                              double new_value,
-                              double old_value);
+                              double new_value, double old_value);
   // Clear a constraint from all its terms.
   virtual void ClearConstraint(MPConstraint* const constraint);
   // Change a coefficient in the linear objective
@@ -166,9 +162,7 @@ class GLPKInterface : public MPSolverInterface {
     return StringPrintf("GLPK %s", glp_version());
   }
 
-  virtual void* underlying_solver() {
-    return reinterpret_cast<void*>(lp_);
-  }
+  virtual void* underlying_solver() { return reinterpret_cast<void*>(lp_); }
 
   virtual double ComputeExactConditionNumber() const;
 
@@ -187,8 +181,7 @@ class GLPKInterface : public MPSolverInterface {
   virtual void SetLpAlgorithm(int value);
 
   void ExtractOldConstraints();
-  void ExtractOneConstraint(MPConstraint* const constraint,
-                            int* const indices,
+  void ExtractOneConstraint(MPConstraint* const constraint, int* const indices,
                             double* const coefs);
   // Transforms basis status from GLPK integer code to MPSolver::BasisStatus.
   MPSolver::BasisStatus TransformGLPKBasisStatus(int glpk_basis_status) const;
@@ -196,16 +189,16 @@ class GLPKInterface : public MPSolverInterface {
   // Computes the L1-norm of the current scaled basis.
   // The L1-norm |A| is defined as max_j sum_i |a_ij|
   // This method is available only for continuous problems.
-  double ComputeScaledBasisL1Norm(
-    int num_rows, int num_cols,
-    double* row_scaling_factor, double* column_scaling_factor) const;
+  double ComputeScaledBasisL1Norm(int num_rows, int num_cols,
+                                  double* row_scaling_factor,
+                                  double* column_scaling_factor) const;
 
   // Computes the L1-norm of the inverse of the current scaled
   // basis.
   // This method is available only for continuous problems.
-  double ComputeInverseScaledBasisL1Norm(
-    int num_rows, int num_cols,
-    double* row_scaling_factor, double* column_scaling_factor) const;
+  double ComputeInverseScaledBasisL1Norm(int num_rows, int num_cols,
+                                         double* row_scaling_factor,
+                                         double* column_scaling_factor) const;
 
   glp_prob* lp_;
   bool mip_;
@@ -281,7 +274,7 @@ void GLPKInterface::SetVariableInteger(int var_index, bool integer) {
   if (mip_) {
     if (var_index != kNoIndex) {
       // Not cached if the variable has been extracted.
-     glp_set_col_kind(lp_, var_index, integer ? GLP_IV : GLP_CV);
+      glp_set_col_kind(lp_, var_index, integer ? GLP_IV : GLP_CV);
     } else {
       sync_status_ = MUST_RELOAD;
     }
@@ -316,8 +309,7 @@ void GLPKInterface::SetConstraintBounds(int index, double lb, double ub) {
 
 void GLPKInterface::SetCoefficient(MPConstraint* const constraint,
                                    const MPVariable* const variable,
-                                   double new_value,
-                                   double old_value) {
+                                   double new_value, double old_value) {
   InvalidateSolutionSynchronization();
   // GLPK does not allow to modify one coefficient at a time, so we
   // extract the whole constraint again, if it has been extracted
@@ -405,15 +397,15 @@ void GLPKInterface::ExtractNewVariables() {
 
 // Extract again existing constraints if they contain new variables.
 void GLPKInterface::ExtractOldConstraints() {
-  int max_constraint_size = solver_->ComputeMaxConstraintSize(
-      0, last_constraint_index_);
+  int max_constraint_size =
+      solver_->ComputeMaxConstraintSize(0, last_constraint_index_);
   // The first entry in the following arrays is dummy, to be
   // consistent with glpk API.
   std::unique_ptr<int[]> indices(new int[max_constraint_size + 1]);
   std::unique_ptr<double[]> coefs(new double[max_constraint_size + 1]);
 
   for (int i = 0; i < last_constraint_index_; ++i) {
-    MPConstraint* const  ct = solver_->constraints_[i];
+    MPConstraint* const ct = solver_->constraints_[i];
     DCHECK_NE(kNoIndex, ct->index());
     const int size = ct->coefficients_.size();
     if (size == 0) {
@@ -441,7 +433,7 @@ void GLPKInterface::ExtractOneConstraint(MPConstraint* const constraint,
     coefs[k] = entry.second;
     ++k;
   }
-  glp_set_mat_row(lp_, constraint->index(), k-1, indices, coefs);
+  glp_set_mat_row(lp_, constraint->index(), k - 1, indices, coefs);
 }
 
 // Define new constraints on old and new variables.
@@ -595,15 +587,13 @@ MPSolver::ResultStatus GLPKInterface::Solve(const MPSolverParameters& param) {
     if (mip_) {
       const double row_activity = glp_mip_row_val(lp_, ct->index());
       ct->set_activity(row_activity);
-      VLOG(4) << "row " << ct->index()
-              << ": activity = " << row_activity;
+      VLOG(4) << "row " << ct->index() << ": activity = " << row_activity;
     } else {
       const double row_activity = glp_get_row_prim(lp_, ct->index());
       ct->set_activity(row_activity);
       const double dual_value = glp_get_row_dual(lp_, ct->index());
       ct->set_dual_value(dual_value);
-      VLOG(4) << "row " << ct->index()
-              << ": activity = " << row_activity
+      VLOG(4) << "row " << ct->index() << ": activity = " << row_activity
               << ": dual value = " << dual_value;
     }
   }
@@ -633,8 +623,7 @@ MPSolver::ResultStatus GLPKInterface::Solve(const MPSolverParameters& param) {
       result_status_ = MPSolver::OPTIMAL;
     } else if (tmp_status == GLP_FEAS) {
       result_status_ = MPSolver::FEASIBLE;
-    } else if (tmp_status == GLP_NOFEAS ||
-               tmp_status == GLP_INFEAS) {
+    } else if (tmp_status == GLP_NOFEAS || tmp_status == GLP_INFEAS) {
       // For infeasible problems, GLPK actually seems to return
       // GLP_UNDEF. So this is never (?) reached.  Return infeasible
       // in case GLPK returns a correct status in future versions.
@@ -654,8 +643,8 @@ MPSolver::ResultStatus GLPKInterface::Solve(const MPSolverParameters& param) {
   return result_status_;
 }
 
-MPSolver::BasisStatus
-GLPKInterface::TransformGLPKBasisStatus(int glpk_basis_status) const {
+MPSolver::BasisStatus GLPKInterface::TransformGLPKBasisStatus(
+    int glpk_basis_status) const {
   switch (glpk_basis_status) {
     case GLP_BS:
       return MPSolver::BASIC;
@@ -676,17 +665,12 @@ GLPKInterface::TransformGLPKBasisStatus(int glpk_basis_status) const {
 // ------ Query statistics on the solution and the solve ------
 
 int64 GLPKInterface::iterations() const {
-  if (mip_) {
-    LOG(WARNING) << "Total number of iterations is not available";
-    return kUnknownNumberOfIterations;
-  } else {
-    if (!CheckSolutionIsSynchronized()) return kUnknownNumberOfIterations;
-#if GLP_MINOR_VERSION >= 49
-    return kUnknownNumberOfIterations;
-#else
+#if GLP_MINOR_VERSION < 49
+  if (!mip_ && CheckSolutionIsSynchronized())
     return lpx_get_int_parm(lp_, LPX_K_ITCNT);
 #endif
-  }
+  LOG(WARNING) << "Total number of iterations is not available";
+  return kUnknownNumberOfIterations;
 }
 
 int64 GLPKInterface::nodes() const {
@@ -732,7 +716,6 @@ MPSolver::BasisStatus GLPKInterface::column_status(int variable_index) const {
   return TransformGLPKBasisStatus(glpk_basis_status);
 }
 
-
 bool GLPKInterface::CheckSolutionExists() const {
   if (result_status_ == MPSolver::ABNORMAL) {
     LOG(WARNING) << "Ignoring ABNORMAL status from GLPK: This status may or may"
@@ -759,9 +742,8 @@ bool GLPKInterface::CheckBestObjectiveBoundExists() const {
 double GLPKInterface::ComputeExactConditionNumber() const {
   if (!IsContinuous()) {
     // TODO(user): support MIP.
-    LOG(DFATAL)
-        << "ComputeExactConditionNumber not implemented for"
-        << " GLPK_MIXED_INTEGER_PROGRAMMING";
+    LOG(DFATAL) << "ComputeExactConditionNumber not implemented for"
+                << " GLPK_MIXED_INTEGER_PROGRAMMING";
     return 0.0;
   }
   if (!CheckSolutionIsSynchronized()) return 0.0;
@@ -779,18 +761,16 @@ double GLPKInterface::ComputeExactConditionNumber() const {
   for (int col = 1; col <= num_cols; ++col) {
     column_scaling_factor[col] = glp_get_sjj(lp_, col);
   }
-  return
-      ComputeInverseScaledBasisL1Norm(
-          num_rows, num_cols,
-          row_scaling_factor.get(), column_scaling_factor.get()) *
-      ComputeScaledBasisL1Norm(
-          num_rows, num_cols,
-          row_scaling_factor.get(), column_scaling_factor.get());
+  return ComputeInverseScaledBasisL1Norm(num_rows, num_cols,
+                                         row_scaling_factor.get(),
+                                         column_scaling_factor.get()) *
+         ComputeScaledBasisL1Norm(num_rows, num_cols, row_scaling_factor.get(),
+                                  column_scaling_factor.get());
 }
 
 double GLPKInterface::ComputeScaledBasisL1Norm(
-    int num_rows, int num_cols,
-    double* row_scaling_factor, double* column_scaling_factor) const {
+    int num_rows, int num_cols, double* row_scaling_factor,
+    double* column_scaling_factor) const {
   double norm = 0.0;
   std::unique_ptr<double[]> values(new double[num_rows + 1]);
   std::unique_ptr<int[]> indices(new int[num_rows + 1]);
@@ -826,8 +806,8 @@ double GLPKInterface::ComputeScaledBasisL1Norm(
 }
 
 double GLPKInterface::ComputeInverseScaledBasisL1Norm(
-    int num_rows, int num_cols,
-    double* row_scaling_factor, double* column_scaling_factor) const {
+    int num_rows, int num_cols, double* row_scaling_factor,
+    double* column_scaling_factor) const {
   // Compute the LU factorization if it doesn't exist yet.
   if (!glp_bf_exists(lp_)) {
     const int factorize_status = glp_factorize(lp_);
@@ -929,7 +909,7 @@ void GLPKInterface::ConfigureGLPKParameters(const MPSolverParameters& param) {
   glp_scale_prob(lp_, GLP_SF_AUTO);
 
   // Use advanced initial basis (options: standard / advanced / Bixby's).
-  glp_adv_basis(lp_, nullptr);
+  glp_adv_basis(lp_, 0);
 
   // Set parameters specified by the user.
   SetParameters(param);
@@ -955,9 +935,7 @@ void GLPKInterface::SetPrimalTolerance(double value) {
   lp_param_.tol_bnd = value;
 }
 
-void GLPKInterface::SetDualTolerance(double value) {
-  lp_param_.tol_dj = value;
-}
+void GLPKInterface::SetDualTolerance(double value) { lp_param_.tol_dj = value; }
 
 void GLPKInterface::SetPresolveMode(int value) {
   switch (value) {
@@ -1004,6 +982,5 @@ MPSolverInterface* BuildGLPKInterface(bool mip, MPSolver* const solver) {
   return new GLPKInterface(solver, mip);
 }
 
-
 }  // namespace operations_research
 #endif  //  #if defined(USE_GLPK)

src/linear_solver/linear_solver.cc

@@ -35,7 +35,9 @@
 #include "base/hash.h"
 #include "base/accurate_sum.h"
 #include "linear_solver/linear_solver.pb.h"
+#include "linear_solver/linear_solver2.pb.h"
 #include "linear_solver/model_exporter.h"
+#include "linear_solver/proto_tools.h"
 #include "util/fp_utils.h"
 
 
@@ -553,6 +555,47 @@ MPSolver::LoadStatus MPSolver::LoadModel(const MPModelProto& input_model) {
   return MPSolver::NO_ERROR;
 }
 
+MPSolver::LoadStatus MPSolver::LoadModelFromProto(
+    const new_proto::MPModelProto& input_model) {
+  for (int i = 0; i < input_model.variable_size(); ++i) {
+    const new_proto::MPVariableProto& var_proto = input_model.variable(i);
+
+    // TODO(user): The variable name var_proto.name() is lost at this point.
+    // This is because MPSolver has no notion of name, just of ids, and these
+    // must be unique which is not necessarily the case of names in the new
+    // proto. MakeNumVar() will just assign an automated variable id when the
+    // the given name argument is empty.
+    //
+    // Note(user): The auto assigned id limits the number of variables to 10^9.
+    MPVariable* variable = MakeNumVar(
+        var_proto.lower_bound(), var_proto.upper_bound(), /*name=*/ "");
+    variable->SetInteger(var_proto.is_integer());
+    SetObjectiveCoefficient(variable, var_proto.objective_coefficient());
+  }
+
+  for (int i = 0; i < input_model.constraint_size(); ++i) {
+    const new_proto::MPConstraintProto& ct_proto = input_model.constraint(i);
+    MPConstraint* const ct = MakeRowConstraint(ct_proto.lower_bound(),
+                                               ct_proto.upper_bound(),
+                                               ct_proto.name());
+    for (int j = 0; j < ct_proto.linear_term_size(); ++j) {
+      const new_proto::MPConstraintProto::UnaryTerm& term_proto =
+          ct_proto.linear_term(j);
+      if (term_proto.var_index() >= variables_.size()) {
+        LOG(ERROR) << "Variable index out of bound in constraint named "
+                    << ct_proto.name() << ".";
+        return MPSolver::UNKNOWN_VARIABLE_ID;
+      }
+      ct->SetCoefficient(variables_[term_proto.var_index()],
+                         term_proto.coefficient());
+    }
+  }
+  SetOptimizationDirection(input_model.maximize());
+  if (input_model.has_objective_offset()) {
+    MutableObjective()->SetOffset(input_model.objective_offset());
+  }
+  return MPSolver::NO_ERROR;
+}
 
 void MPSolver::FillSolutionResponse(MPSolutionResponse* response) const {
   CHECK_NOTNULL(response);
@@ -612,6 +655,40 @@ void MPSolver::FillSolutionResponse(MPSolutionResponse* response) const {
   }
 }
 
+namespace {
+new_proto::MPSolutionResponse::Status ResultStatusToMPSolutionResponse(
+    MPSolver::ResultStatus status) {
+  switch (status) {
+    case MPSolver::OPTIMAL :
+      return new_proto::MPSolutionResponse::OPTIMAL;
+    case MPSolver::FEASIBLE :
+      return new_proto::MPSolutionResponse::FEASIBLE;
+    case MPSolver::INFEASIBLE :
+      return new_proto::MPSolutionResponse::INFEASIBLE;
+    case MPSolver::UNBOUNDED :
+      return new_proto::MPSolutionResponse::UNBOUNDED;
+    case MPSolver::ABNORMAL :
+      return new_proto::MPSolutionResponse::ABNORMAL;
+    default:
+      return new_proto::MPSolutionResponse::UNKNOWN;
+  }
+}
+}  // namespace
+
+void MPSolver::FillSolutionResponseProto(
+    new_proto::MPSolutionResponse* response) const {
+  CHECK_NOTNULL(response);
+  response->Clear();
+  response->set_status(
+      ResultStatusToMPSolutionResponse(interface_->result_status_));
+  if (interface_->result_status_ == MPSolver::OPTIMAL ||
+      interface_->result_status_ == MPSolver::FEASIBLE) {
+    response->set_objective_value(objective_value());
+    for (int i = 0; i < variables_.size(); ++i) {
+      response->add_variable_value(variables_[i]->solution_value());
+    }
+  }
+}
 
 // static
 void MPSolver::SolveWithProtocolBuffers(const MPModelRequest& model_request,
@@ -624,6 +701,8 @@ void MPSolver::SolveWithProtocolBuffers(const MPModelRequest& model_request,
   const MPSolver::LoadStatus loadStatus = solver.LoadModel(model);
   if (loadStatus != MPSolver::NO_ERROR) {
     LOG(WARNING) << "Loading model from protocol buffer failed, "
+                 << "load status = "
+                 << Error::Code_Name(static_cast<Error::Code>(loadStatus))
                  << " (" << loadStatus << ")";
     response->set_result_status(MPSolutionResponse::ABNORMAL);
     return;
@@ -636,6 +715,33 @@ void MPSolver::SolveWithProtocolBuffers(const MPModelRequest& model_request,
   solver.FillSolutionResponse(response);
 }
 
+// static
+void MPSolver::SolveWithProto(
+    const new_proto::MPModelRequest& model_request,
+    new_proto::MPSolutionResponse* response) {
+  CHECK_NOTNULL(response);
+  const new_proto::MPModelProto& model = model_request.model();
+  MPSolver solver(model.name(),
+                  static_cast<MPSolver::OptimizationProblemType>(
+                      model_request.solver_type()));
+  const MPSolver::LoadStatus loadStatus = solver.LoadModelFromProto(model);
+  if (loadStatus != MPSolver::NO_ERROR) {
+    LOG(WARNING) << "Loading model from protocol buffer failed, "
+                 << "load status = "
+                 << Error::Code_Name(static_cast<Error::Code>(loadStatus))
+                 << " (" << loadStatus << ")";
+    response->set_status(new_proto::MPSolutionResponse::ABNORMAL);
+    return;
+  }
+  if (model_request.has_solver_time_limit_seconds()) {
+    // static_cast<int64> avoids a warning with -Wreal-conversion. This
+    // helps catching bugs with unwanted conversions from double to ints.
+    solver.set_time_limit(
+        static_cast<int64>(model_request.solver_time_limit_seconds()) * 1000);
+  }
+  solver.Solve();
+  solver.FillSolutionResponseProto(response);
+}
 
 namespace {
 // Outputs the terms in var_coeff_map to output_term, by sorting the
@@ -673,32 +779,32 @@ void OutputTermsToProto(
 
 void MPSolver::ExportModel(MPModelProto* output_model) const {
   DCHECK(output_model != NULL);
-  if (output_model->variables_size() > 0 ||
-      output_model->has_maximize() ||
-      output_model->objective_terms_size() > 0 ||
-      output_model->constraints_size() > 0 ||
-      output_model->has_name() ||
-      output_model->has_objective_offset()) {
-    LOG(WARNING) << "The model protocol buffer is not empty, "
-                 << "it will be overwritten.";
-    output_model->clear_variables();
-    output_model->clear_maximize();
-    output_model->clear_objective_terms();
-    output_model->clear_constraints();
-    output_model->clear_name();
-  }
+  new_proto::MPModelProto new_proto_model;
+  ExportModelToNewProto(&new_proto_model);
+  ConvertNewMPModelProtoToOld(new_proto_model, output_model);
+}
 
+void MPSolver::ExportModelToNewProto(
+    new_proto::MPModelProto* output_model) const {
+  DCHECK(output_model != NULL);
+  output_model->Clear();
   // Name
   output_model->set_name(Name());
   // Variables
   for (int j = 0; j < variables_.size(); ++j) {
     const MPVariable* const var = variables_[j];
-    MPVariableProto* const variable_proto = output_model->add_variables();
-    DCHECK(!var->name().empty());
-    variable_proto->set_id(var->name());
-    variable_proto->set_lb(var->lb());
-    variable_proto->set_ub(var->ub());
-    variable_proto->set_integer(var->integer());
+    new_proto::MPVariableProto* const variable_proto
+        = output_model->add_variable();
+    // TODO(user): Add option to avoid filling the var name to avoid overly
+    // large protocol buffers.
+    variable_proto->set_name(var->name());
+    variable_proto->set_lower_bound(var->lb());
+    variable_proto->set_upper_bound(var->ub());
+    variable_proto->set_is_integer(var->integer());
+    if (objective_->GetCoefficient(var) != 0.0) {
+      variable_proto->set_objective_coefficient(
+          objective_->GetCoefficient(var));
+    }
   }
 
   // Map the variables to their indices. This is needed to output the
@@ -707,37 +813,56 @@ void MPSolver::ExportModel(MPModelProto* output_model) const {
   // This step is needed as long as the variable indices are given by the
   // underlying solver at the time of model extraction.
   // TODO(user): remove this step.
-  hash_map<const MPVariable*, int> var_name_to_index;
+  hash_map<const MPVariable*, int> var_to_index;
   for (int j = 0; j < variables_.size(); ++j) {
-    var_name_to_index[variables_[j]] = j;
+    var_to_index[variables_[j]] = j;
   }
 
   // Constraints
   for (int i = 0; i < constraints_.size(); ++i) {
     MPConstraint* const constraint = constraints_[i];
-    MPConstraintProto* const constraint_proto = output_model->add_constraints();
-    // Constraint names need to be non-empty.
-    DCHECK(!constraint->name().empty());
-    constraint_proto->set_id(constraint->name());
-    constraint_proto->set_lb(constraint->lb());
-    constraint_proto->set_ub(constraint->ub());
-    OutputTermsToProto(variables_, var_name_to_index, constraint->coefficients_,
-                       constraint_proto->mutable_terms());
+    new_proto::MPConstraintProto* const constraint_proto
+        = output_model->add_constraint();
+    constraint_proto->set_name(constraint->name());
+    constraint_proto->set_lower_bound(constraint->lb());
+    constraint_proto->set_upper_bound(constraint->ub());
+    constraint_proto->set_is_lazy(constraint->is_lazy());
+    // Vector linear_term will contain pairs (variable index, coeff), that will
+    // be sorted by variable index.
+    std::vector<std::pair<int, double> > linear_term;
+    for (CoeffEntry entry : constraint->coefficients_) {
+      const MPVariable* const var = entry.first;
+      const int var_index = FindWithDefault(var_to_index, var, -1);
+      DCHECK_NE(-1, var_index);
+      const double coeff = entry.second;
+      linear_term.push_back(std::pair<int, double>(var_index, coeff));
+    }
+    // The cost of std::sort is expected to be low as constraints usually have very
+    // few terms.
+    std::sort(linear_term.begin(), linear_term.end());
+    // Now use linear term.
+    for (int k = 0; k < linear_term.size(); ++k) {
+      const std::pair<int, double>& p = linear_term[k];
+      const int var_index = p.first;
+      const double coeff = p.second;
+      new_proto::MPConstraintProto_UnaryTerm* unary_term
+          = constraint_proto->add_linear_term();
+      unary_term->set_var_index(var_index);
+      unary_term->set_coefficient(coeff);
+      constraint_proto->add_var_index(var_index);
+      constraint_proto->add_coefficient(coeff);
+    }
   }
-
-  // Objective.
-  OutputTermsToProto(variables_, var_name_to_index, objective_->coefficients_,
-                     output_model->mutable_objective_terms());
   output_model->set_maximize(Objective().maximization());
   output_model->set_objective_offset(Objective().offset());
 }
 
-
-bool MPSolver::LoadSolutionFromProto(const MPSolutionResponse& response) {
+bool MPSolver::LoadSolutionFromNewProto(
+    const new_proto::MPSolutionResponse& response) {
   interface_->result_status_ =
-      static_cast<ResultStatus>(response.result_status());
-  if (response.result_status() != MPSolutionResponse::OPTIMAL &&
-      response.result_status() != MPSolutionResponse::FEASIBLE) {
+      static_cast<ResultStatus>(response.status());
+  if (response.status() != new_proto::MPSolutionResponse::OPTIMAL &&
+      response.status() != new_proto::MPSolutionResponse::FEASIBLE) {
     LOG(ERROR)
         << "Cannot load a solution unless its status is OPTIMAL or FEASIBLE.";
     return false;
@@ -745,55 +870,37 @@ bool MPSolver::LoadSolutionFromProto(const MPSolutionResponse& response) {
   // Before touching the variables, verify that the solution looks legit:
   // each variable of the MPSolver must have its value listed exactly once, and
   // each listed solution should correspond to a known variable.
-  if (response.solution_values_size() != variables_.size()) {
+  if (response.variable_value_size() != variables_.size()) {
     LOG(ERROR)
         << "Trying to load a solution whose number of variables does not"
         << " correspond to the Solver.";
     return false;
   }
-  std::vector<bool> variable_has_solution_value(variables_.size(), false);
+  double largest_error = 0;
+  interface_->ExtractModel();
   int num_vars_out_of_bounds = 0;
-  for (int i = 0; i < response.solution_values_size(); ++i) {
-    const string& var_name = response.solution_values(i).variable_id();
-    const int var_index =
-        FindWithDefault(variable_name_to_index_, var_name, -1);
-    if (var_index < 0) {
-      LOG(ERROR)
-          << "Trying to load a solution with unknown var '" << var_name << "'.";
-      return false;
-    }
-    if (variable_has_solution_value[var_index]) {
-      LOG(ERROR)
-          << "Trying to load a solution value where variable '" << var_name
-          << "' has its solution value listed twice.";
-      return false;
-    }
-    variable_has_solution_value[var_index] = true;
+  const double tolerance = MPSolverParameters::kDefaultPrimalTolerance;
+  for (int i = 0; i< response.variable_value_size(); ++i) {
     // Look further: verify the bounds. Since linear solvers yield (small)
     // numerical errors, though, we just log a warning if the variables look
     // like they are out of their bounds. The user should inspect the values.
-    const double var_value = response.solution_values(i).value();
-    const MPVariable* var = variables_[var_index];
-    DCHECK_EQ(var->name(), var_name);
-    if (var_value < var->lb() || var_value > var->ub()) {
+    const double var_value = response.variable_value(i);
+    MPVariable* var = variables_[i];
+    // TODO(user): Use parameter when they become available in this class.
+    const double lb_error = var->lb() - var_value;
+    const double ub_error = var_value - var->ub();
+    if (lb_error > tolerance || ub_error > tolerance) {
       ++num_vars_out_of_bounds;
+      largest_error = std::max(largest_error, std::max(lb_error, ub_error));
     }
+    var->set_solution_value(var_value);
   }
   if (num_vars_out_of_bounds > 0) {
     LOG(WARNING)
         << "Loaded a solution whose variables matched the solver's, but "
         << num_vars_out_of_bounds << " out of " << variables_.size()
-        << " had values outside their bounds.";
-  }
-  // At this point, we know for sure that the solution can be safely loaded.
-  // We extract the model onto the underlying solver, so that the accessors work
-  // properly later (as if we had called Solve()).
-  interface_->ExtractModel();
-  // Do a second pass, where we actually set the variable values.
-  for (int i = 0; i < response.solution_values_size(); ++i) {
-    variables_[FindOrDie(variable_name_to_index_,
-                         response.solution_values(i).variable_id())]
-        ->set_solution_value(response.solution_values(i).value());
+        << " exceed one of their bounds by more tahn the primal tolerance: "
+        << tolerance;
   }
   // Set the objective value, if is known.
   if (response.has_objective_value()) {
@@ -1579,5 +1686,81 @@ int MPSolverParameters::GetIntegerParam(
 }
 
 
+bool MPSolver::LoadSolutionFromProto(const MPSolutionResponse& response) {
+  // CHANGES TO THIS METHOD ARE DISCOURAGED.
+  // This method will be deprecated in 2014.
+  // Please try using LoadSolutionFromNewProto() instead.
+  // TODO(user): Deprecate then remove this code.
+  interface_->result_status_ =
+      static_cast<ResultStatus>(response.result_status());
+  if (response.result_status() != MPSolutionResponse::OPTIMAL &&
+      response.result_status() != MPSolutionResponse::FEASIBLE) {
+    LOG(ERROR)
+        << "Cannot load a solution unless its status is OPTIMAL or FEASIBLE.";
+    return false;
+  }
+  // Before touching the variables, verify that the solution looks legit:
+  // each variable of the MPSolver must have its value listed exactly once, and
+  // each listed solution should correspond to a known variable.
+  if (response.solution_values_size() != variables_.size()) {
+    LOG(ERROR)
+        << "Trying to load a solution whose number of variables does not"
+        << " correspond to the Solver.";
+    return false;
+  }
+  std::vector<bool> variable_has_solution_value(variables_.size(), false);
+  int num_vars_out_of_bounds = 0;
+  for (int i = 0; i < response.solution_values_size(); ++i) {
+    const string& var_name = response.solution_values(i).variable_id();
+    const int var_index =
+        FindWithDefault(variable_name_to_index_, var_name, -1);
+    if (var_index < 0) {
+      LOG(ERROR)
+          << "Trying to load a solution with unknown var '" << var_name << "'.";
+      return false;
+    }
+    if (variable_has_solution_value[var_index]) {
+      LOG(ERROR)
+          << "Trying to load a solution value where variable '" << var_name
+          << "' has its solution value listed twice.";
+      return false;
+    }
+    variable_has_solution_value[var_index] = true;
+    // Look further: verify the bounds. Since linear solvers yield (small)
+    // numerical errors, though, we just log a warning if the variables look
+    // like they are out of their bounds. The user should inspect the values.
+    const double var_value = response.solution_values(i).value();
+    const MPVariable* var = variables_[var_index];
+    DCHECK_EQ(var->name(), var_name);
+    if (var_value < var->lb() || var_value > var->ub()) {
+      ++num_vars_out_of_bounds;
+    }
+  }
+  if (num_vars_out_of_bounds > 0) {
+    LOG(WARNING)
+        << "Loaded a solution whose variables matched the solver's, but "
+        << num_vars_out_of_bounds << " out of " << variables_.size()
+        << " had values outside their bounds.";
+  }
+  // At this point, we know for sure that the solution can be safely loaded.
+  // We extract the model onto the underlying solver, so that the accessors work
+  // properly later (as if we had called Solve()).
+  interface_->ExtractModel();
+  // Do a second pass, where we actually set the variable values.
+  for (int i = 0; i < response.solution_values_size(); ++i) {
+    variables_[FindOrDie(variable_name_to_index_,
+                         response.solution_values(i).variable_id())]
+        ->set_solution_value(response.solution_values(i).value());
+  }
+  // Set the objective value, if is known.
+  if (response.has_objective_value()) {
+    interface_->objective_value_ = response.objective_value();
+  }
+  // Mark the status as SOLUTION_SYNCHRONIZED, so that users may inspect the
+  // solution normally.
+  interface_->sync_status_ = MPSolverInterface::SOLUTION_SYNCHRONIZED;
+  return true;
+}
+
 }  // namespace operations_research
 

src/linear_solver/linear_solver.h

@@ -158,6 +158,13 @@ class MPSolverInterface;
 class MPSolverParameters;
 class MPVariable;
 
+// The new MP protocol buffer format. New clients that want to work with
+// protocol buffers should use this.
+namespace new_proto {
+  class MPModelProto;
+  class MPModelRequest;
+  class MPSolutionResponse;
+}  // namespace new_proto
 
 // This mathematical programming (MP) solver class is the main class
 // though which users build and solve problems.
@@ -351,12 +358,12 @@ class MPSolver {
   };
 
   // Loads model from protocol buffer.
-  LoadStatus LoadModel(const MPModelProto& input_model);
+  LoadStatus LoadModelFromProto(const new_proto::MPModelProto& input_model);
 
   // Encodes the current solution in a solution response protocol buffer.
   // Only nonzero variable values are stored in order to reduce the
   // size of the MPSolutionResponse protocol buffer.
-  void FillSolutionResponse(MPSolutionResponse* response) const;
+  void FillSolutionResponseProto(new_proto::MPSolutionResponse* response) const;
 
   // Solves the model encoded by a MPModelRequest protocol buffer and
   // fills the solution encoded as a MPSolutionResponse.
@@ -364,27 +371,28 @@ class MPSolver {
   // end. If you want to keep the MPSolver alive (for debugging, or for
   // incremental solving), you should write another version of this function
   // that creates the MPSolver object on the heap and returns it.
-  static void SolveWithProtocolBuffers(const MPModelRequest& model_request,
-                                       MPSolutionResponse* response);
-
+  static void SolveWithProto(
+      const new_proto::MPModelRequest& model_request,
+      new_proto::MPSolutionResponse* response);
 
   // Exports model to protocol buffer.
-  void ExportModel(MPModelProto* output_model) const;
+  // TODO(user): rename to ExportModelToProto when possible.
+  void ExportModelToNewProto(new_proto::MPModelProto* output_model) const;
 
-
-  // Load a solution encoded in a protocol buffer onto this solver.
+  // Load a solution encoded in a protocol buffer onto this solver for easy
+  // access via the MPSolver interface.
   //
   // IMPORTANT: This may only be used in conjunction with ExportModel(),
   // following this example:
   //   MPSolver my_solver;
   //   ... add variables and constraints ...
-  //   MPModelProto model_proto;
-  //   my_solver.ExportModel(&model_proto);
-  //   MPSolutionResponse solver_response;
+  //   new_proto::MPModelProto model_proto;
+  //   my_solver.ExportModelToNewProto(&model_proto);
+  //   new_proto::MPSolutionResponse solver_response;
   //   // This can be replaced by a stubby call to the linear solver server.
-  //   MPSolver::SolveWithProtocolBuffers(model_proto, &solver_response);
+  //   MPSolver::SolveWithProto(model_proto, &solver_response);
   //   if (solver_response.result_status() == MPSolutionResponse::OPTIMAL) {
-  //     CHECK(my_solver.LoadSolutionFromProto(solver_response));
+  //     CHECK(my_solver.LoadSolutionFromNewProto(solver_response));
   //     ... inspect the solution using the usual API: solution_value(), etc...
   //   }
   //
@@ -400,12 +408,14 @@ class MPSolver {
   // - loading a solution with a status other than OPTIMAL / FEASIBLE.
   // Note: the variable and objective values aren't checked. You can use
   // VerifySolution() for that.
-  bool LoadSolutionFromProto(const MPSolutionResponse& response);
+  // TODO(user): rename to LoadSolutionFromProto when possible.
+  bool LoadSolutionFromNewProto(
+      const new_proto::MPSolutionResponse& response);
 
   // ----- Export model to files or strings -----
 
   // Shortcuts to the homonymous MPModelProtoExporter methods, via
-  // exporting to a MPModelProto with ExportModel() (see above).
+  // exporting to a MPModelProto with ExportModelToNewProto() (see above).
   bool ExportModelAsLpFormat(bool obfuscated, string* model_str);
   bool ExportModelAsMpsFormat(
       bool fixed_format, bool obfuscated, string* model_str);
@@ -548,6 +558,16 @@ class MPSolver {
   bool Maximization() const;
   bool Minimization() const;
 
+  // DEPRECATED methods. Use the explicitly listed replacement.
+  LoadStatus LoadModel(const MPModelProto& input_model);  // LoadModelFromProto
+  void FillSolutionResponse(  // FillSolutionResponseProto
+      MPSolutionResponse* response) const;
+  static void SolveWithProtocolBuffers(  // SolveWithProto
+      const MPModelRequest& model_request, MPSolutionResponse* response);
+  void ExportModel(MPModelProto* output_model) const;  // ExportModelToNewProto
+  bool LoadSolutionFromProto(  // LoadSolutionFromNewProto
+      const MPSolutionResponse& response);
+
  private:
   // Computes the size of the constraint with the largest number of
   // coefficients with index in [min_constraint_index,

src/linear_solver/linear_solver.proto

@@ -87,3 +87,14 @@ message MPSolutionResponse {
   repeated MPSolutionValue solution_values = 3;
 }
 
+message Error {
+  enum Code {
+    NO_ERROR = 0;
+    INVALID_PROBLEM_TYPE = 1;
+    DUPLICATE_VARIABLE_ID = 2;
+    UNKNOWN_VARIABLE_ID = 3;
+    REQUEST_IS_QOD = 4;
+    RPC_DEADLINE_TOO_SMALL = 5;
+  }
+}
+

src/linear_solver/linear_solver.swig

@@ -20,6 +20,7 @@
 %{
 #include "linear_solver/linear_solver.h"
 #include "linear_solver/linear_solver.pb.h"
+#include "linear_solver/linear_solver2.pb.h"
 #include "linear_solver/linear_solver_ext.h"
 %}
 
@@ -526,6 +527,9 @@ namespace operations_research {
 // The following 3 methods that use protocol buffers as output
 // arguments are replaced by methods that return a protocol buffer,
 // see code below.
+%ignore MPSolver::ExportModelToNewProto;
+%ignore MPSolver::FillSolutionResponseProto;
+%ignore MPSolver::SolveWithProto;
 %ignore MPSolver::ExportModel;
 %ignore MPSolver::FillSolutionResponse;
 %ignore MPSolver::SolveWithProtocolBuffers;
@@ -635,6 +639,7 @@ namespace operations_research {
 // arguments are replaced by methods that return a protocol buffer,
 // see code below.
 %ignore MPSolver::ExportModel;
+%ignore MPSolver::ExportModelToNewProto;
 %ignore MPSolver::FillSolutionResponse;
 %ignore MPSolver::SolveWithProtocolBuffers;