more pdlp work

ortools/pdlp/primal_dual_hybrid_gradient.cc

@@ -144,8 +144,9 @@ std::string ConvergenceInformationString(
           convergence_information.l2_primal_variable(),
           convergence_information.l2_dual_variable());
     case OPTIMALITY_NORM_UNSPECIFIED:
-      LOG(FATAL) << "Invalid residual norm.";
+      LOG(FATAL) << "Residual norm not specified.";
   }
+  LOG(FATAL) << "Invalid residual norm " << residual_norm << ".";
 }
 
 std::string ConvergenceInformationShortString(
@@ -178,8 +179,9 @@ std::string ConvergenceInformationShortString(
           convergence_information.primal_objective(),
           convergence_information.dual_objective());
     case OPTIMALITY_NORM_UNSPECIFIED:
-      LOG(FATAL) << "Invalid residual norm.";
+      LOG(FATAL) << "Residual norm not specified.";
   }
+  LOG(FATAL) << "Invalid residual norm " << residual_norm << ".";
 }
 
 void LogInfoWithoutPrefix(absl::string_view message) {

ortools/pdlp/quadratic_program_io.cc

@@ -29,16 +29,18 @@
 #include "absl/status/status.h"
 #include "absl/status/statusor.h"
 #include "absl/strings/match.h"
+#include "absl/strings/str_format.h"
 #include "ortools/base/basictypes.h"
 #include "ortools/base/file.h"
 #include "ortools/base/helpers.h"
 #include "ortools/base/logging.h"
 #include "ortools/base/mathutil.h"
 #include "ortools/base/options.h"
+#include "ortools/base/status_builder.h"
 #include "ortools/base/status_macros.h"
 #include "ortools/linear_solver/linear_solver.pb.h"
 #include "ortools/linear_solver/model_exporter.h"
-#include "ortools/lp_data/mps_reader.h"
+#include "ortools/lp_data/mps_reader_template.h"
 #include "ortools/pdlp/quadratic_program.h"
 #include "ortools/util/file_util.h"
 
@@ -64,21 +66,6 @@ QuadraticProgram ReadQuadraticProgramOrDie(const std::string& filename,
               << ".json, and .json.gz";
 }
 
-QuadraticProgram ReadMpsLinearProgramOrDie(const std::string& lp_file,
-                                           bool include_names) {
-  absl::StatusOr<MPModelProto> lp_proto =
-      operations_research::glop::MpsFileToMPModelProto(lp_file);
-  QCHECK_OK(lp_proto);
-  // `MpsFileToMPModelProto` sometimes fails silently if the file isn't read
-  // properly.
-  QCHECK_GT(lp_proto->variable_size(), 0)
-      << "No variables in LP. Error reading file? " << lp_file;
-  auto result = QpFromMpModelProto(*lp_proto, /*relax_integer_variables=*/true,
-                                   include_names);
-  QCHECK_OK(result);
-  return *std::move(result);
-}
-
 QuadraticProgram ReadMPModelProtoFileOrDie(
     const std::string& mpmodel_proto_file, bool include_names) {
   MPModelProto lp_proto;
@@ -112,4 +99,311 @@ absl::Status WriteQuadraticProgramToMPModelProto(
   return file::SetBinaryProto(mpmodel_proto_file, proto, file::Defaults());
 }
 
+namespace {
+
+// Class implementing the
+// `ortools/lp_data/mps_reader_template.h` interface that only
+// stores the names of rows and columns, and the number of non-zeros found.
+class MpsReaderDimensionAndNames {
+ public:
+  using IndexType = int64_t;
+  void SetUp() {
+    read_or_parse_failed_ = false;
+    col_name_to_index_.clear();
+    row_name_to_index_.clear();
+    added_non_zeros_ = 0;
+  }
+  void CleanUp() {}
+  double ConstraintLowerBound(IndexType index) { return 0; }
+  double ConstraintUpperBound(IndexType index) { return 0; }
+  IndexType FindOrCreateConstraint(const std::string& row_name) {
+    const auto [iterator, success_unused] =
+        row_name_to_index_.insert({row_name, row_name_to_index_.size()});
+    return iterator->second;
+  }
+  IndexType FindOrCreateVariable(const std::string& col_name) {
+    const auto [iterator, success_unused] =
+        col_name_to_index_.insert({col_name, col_name_to_index_.size()});
+    return iterator->second;
+  }
+  void SetConstraintBounds(IndexType row_index, double lower_bound,
+                           double upper_bound) {}
+  void SetConstraintCoefficient(IndexType row_index, IndexType col_index,
+                                double coefficient) {
+    ++added_non_zeros_;
+  }
+  void SetIsLazy(IndexType row_index) {}
+  void SetName(const std::string& problem_name) {}
+  void SetObjectiveCoefficient(IndexType index, double coefficient) {}
+  void SetObjectiveDirection(bool maximize) {}
+  void SetObjectiveOffset(double offset) {}
+  void SetVariableTypeToInteger(IndexType index) {}
+  void SetVariableTypeToSemiContinuous(IndexType index) {
+    LOG_FIRST_N(WARNING, 1)
+        << "Semi-continuous variables not supported, failed to parse file";
+    read_or_parse_failed_ = true;
+  }
+  void SetVariableBounds(IndexType index, double lower_bound,
+                         double upper_bound) {}
+  double VariableLowerBound(IndexType index) { return 0; }
+  double VariableUpperBound(IndexType index) { return 0; }
+  absl::Status CreateIndicatorConstraint(const std::string& row_name,
+                                         IndexType col_index, bool var_value) {
+    absl::string_view message =
+        "Indicator constraints not supported, failed to parse file";
+    LOG_FIRST_N(WARNING, 1) << message;
+    return absl::InvalidArgumentError(message);
+  }
+
+  // Lookup constant functions. They assume that `row_name` or `col_name`,
+  // respectively, have been previously added to the internal set of rows or
+  // variables. The functions QFAIL if this is not true.
+  IndexType FindVariableOrDie(const std::string& col_name) const {
+    const auto it = col_name_to_index_.find(col_name);
+    if (it == col_name_to_index_.end()) {
+      LOG(QFATAL) << absl::StrFormat(
+          "column `%s` not previously found in file, internal error?",
+          col_name);
+    }
+    return it->second;
+  }
+  IndexType FindConstraintOrDie(const std::string& row_name) const {
+    const auto it = row_name_to_index_.find(row_name);
+    if (it == row_name_to_index_.end()) {
+      LOG(QFATAL) << absl::StrFormat(
+          "row `%s` not previously found in file, internal error?", row_name);
+    }
+    return it->second;
+  }
+
+  // Number of non-zeros added so-far.
+  int64_t AddedNonZeros() const { return added_non_zeros_; }
+
+  // Number of variables added so-far.
+  int64_t NumVariables() const { return col_name_to_index_.size(); }
+
+  // Number of constraints added so-far.
+  int64_t NumConstraints() const { return row_name_to_index_.size(); }
+
+  // Return `true` if an unsupported MPS section was found.
+  bool FailedToParse() const { return read_or_parse_failed_; }
+
+  // Return a const-hash map of col names to indices.
+  const absl::flat_hash_map<std::string, IndexType>& ColNameIndexMap() const {
+    return col_name_to_index_;
+  }
+
+  // Return a const-hash map of row names to indices.
+  const absl::flat_hash_map<std::string, IndexType>& RowNameIndexMap() const {
+    return row_name_to_index_;
+  }
+
+ private:
+  bool read_or_parse_failed_ = false;
+  absl::flat_hash_map<std::string, IndexType> col_name_to_index_;
+  absl::flat_hash_map<std::string, IndexType> row_name_to_index_;
+  int64_t added_non_zeros_ = 0;
+};
+
+// Class implementing the
+// `ortools/lp_data/mps_reader_template.h` interface.
+// It is intended to be used in conjunction with `MpsReaderDimensionAndNames`
+// as follows:
+//
+// // Retrieve dimension and name information from file.
+// MpsReaderDimensionAndNames dimension_and_name;
+// MPSReaderTemplate<MpsReaderDimensionAndNames> dimension_parser;
+// dimension_parser.ParseFile(file_name, &dimension_and_name);
+// // Store QP problem coefficients.
+// MpsReaderQpDataWrapper qp_wrapper(&dimension_and_name, include_names);
+// MPSReaderTemplate<MpsReaderQpDataWrapper> qp_parser;
+// qp_parser.ParseFile(file_name, &qp_wrapper);
+// // Retrieve fully assembled QP.
+// QuadraticProgram result = qp_wrapper.GetAndClearQuadraticProgram();
+class MpsReaderQpDataWrapper {
+ public:
+  using IndexType = int64_t;
+  // Constructor, `*dimension_and_names` must outlive the class, and be
+  // constant during the object's lifetime. If `include_names`=true, the
+  // resulting QuadraticProgram from `GetAndClearQuadraticProgram()` will
+  // include name information from the MPS file.
+  // NOTE: The code assumes that the file to be read is the same file already
+  // read using the `*dimension_and_names` argument.
+  MpsReaderQpDataWrapper(const MpsReaderDimensionAndNames* dimension_and_names,
+                         bool include_names)
+      : include_names_{include_names},
+        dimension_and_names_{*dimension_and_names} {}
+
+  // Implements the `MPSReaderTemplate` interface.
+  void SetUp() {
+    const int64_t num_variables = dimension_and_names_.NumVariables();
+    const int64_t num_constraints = dimension_and_names_.NumConstraints();
+    triplets_.reserve(dimension_and_names_.AddedNonZeros());
+    quadratic_program_ = QuadraticProgram(/*num_variables=*/num_variables,
+                                          /*num_constraints=*/num_constraints);
+    // Default variables in MPS files have a zero lower bound, an infinity
+    // upper bound, and a zero objective; while default constraints are
+    // 'equal to zero' constraints.
+    quadratic_program_.constraint_lower_bounds =
+        Eigen::VectorXd::Zero(num_constraints);
+    quadratic_program_.constraint_upper_bounds =
+        Eigen::VectorXd::Zero(num_constraints);
+    quadratic_program_.variable_lower_bounds =
+        Eigen::VectorXd::Zero(num_variables);
+  }
+  void CleanUp() {
+    SetEigenMatrixFromTriplets(std::move(triplets_),
+                               quadratic_program_.constraint_matrix);
+    // Deal with maximization problems.
+    if (quadratic_program_.objective_scaling_factor == -1) {
+      quadratic_program_.objective_offset *= -1;
+      for (double& objective_coefficient :
+           quadratic_program_.objective_vector) {
+        objective_coefficient *= -1;
+      }
+    }
+    if (include_names_) {
+      quadratic_program_.variable_names =
+          std::vector<std::string>(dimension_and_names_.NumVariables());
+      quadratic_program_.constraint_names =
+          std::vector<std::string>(dimension_and_names_.NumConstraints());
+      for (auto [name, index] : dimension_and_names_.ColNameIndexMap()) {
+        (*quadratic_program_.variable_names)[index] = name;
+      }
+      for (auto [name, index] : dimension_and_names_.RowNameIndexMap()) {
+        (*quadratic_program_.constraint_names)[index] = name;
+      }
+    }
+  }
+  double ConstraintLowerBound(IndexType index) {
+    return quadratic_program_.constraint_lower_bounds[index];
+  }
+  double ConstraintUpperBound(IndexType index) {
+    return quadratic_program_.constraint_upper_bounds[index];
+  }
+  IndexType FindOrCreateConstraint(const std::string& row_name) {
+    return dimension_and_names_.FindConstraintOrDie(row_name);
+  }
+  IndexType FindOrCreateVariable(const std::string& col_name) {
+    return dimension_and_names_.FindVariableOrDie(col_name);
+  }
+  void SetConstraintBounds(IndexType index, double lower_bound,
+                           double upper_bound) {
+    quadratic_program_.constraint_lower_bounds[index] = lower_bound;
+    quadratic_program_.constraint_upper_bounds[index] = upper_bound;
+  }
+  void SetConstraintCoefficient(IndexType row_index, IndexType col_index,
+                                double coefficient) {
+    DCHECK_LT(triplets_.size(), triplets_.capacity());
+    triplets_.emplace_back(row_index, col_index, coefficient);
+  }
+  void SetIsLazy(IndexType row_index) {
+    LOG_FIRST_N(WARNING, 1) << "Lazy constraint information lost, treated as "
+                               "regular constraint instead";
+  }
+  void SetName(const std::string& problem_name) {
+    if (include_names_) {
+      quadratic_program_.problem_name = problem_name;
+    }
+  }
+  void SetObjectiveCoefficient(IndexType index, double coefficient) {
+    quadratic_program_.objective_vector[index] = coefficient;
+  }
+  void SetObjectiveDirection(bool maximize) {
+    quadratic_program_.objective_scaling_factor = (maximize ? -1 : 1);
+  }
+  void SetObjectiveOffset(double offset) {
+    quadratic_program_.objective_offset = offset;
+  }
+  void SetVariableTypeToInteger(IndexType index) {
+    LOG_FIRST_N(WARNING, 1) << "Dropping integrality requirements, all "
+                               "variables treated as continuous";
+  }
+  void SetVariableTypeToSemiContinuous(IndexType index) {
+    // This line should never execute, as the code must fail on
+    // `MpsReaderDimensionAndNames::SetVariableTypeToSemiContinuous` before.
+    LOG(QFATAL) << "Semi-continuous variables are unsupported, and should have "
+                   "been detected before (in MpsReaderDimensionAndNames)";
+  }
+  void SetVariableBounds(IndexType index, double lower_bound,
+                         double upper_bound) {
+    quadratic_program_.variable_lower_bounds[index] = lower_bound;
+    quadratic_program_.variable_upper_bounds[index] = upper_bound;
+  }
+  double VariableLowerBound(IndexType index) {
+    return quadratic_program_.variable_lower_bounds[index];
+  }
+  double VariableUpperBound(IndexType index) {
+    return quadratic_program_.variable_upper_bounds[index];
+  }
+  absl::Status CreateIndicatorConstraint(const std::string& row_name,
+                                         IndexType col_index, bool var_value) {
+    // This function should never be called, as the code must fail on
+    // `MpsReaderDimensionAndNames::CreateIndicatorConstraint` before.
+    LOG(QFATAL) << "Indicator constraints are unsupported, and should have "
+                   "been detected before (in MpsReaderDimensionAndNames)";
+  }
+
+  // Returns a `QuadraticProgram` holding all information read by the
+  // `mps_reader_template.h` interface. It leaves the internal quadratic
+  // program in an indeterminate state.
+  QuadraticProgram&& GetAndClearQuadraticProgram() {
+    return std::move(quadratic_program_);
+  }
+
+ private:
+  bool include_names_;
+  QuadraticProgram quadratic_program_;
+  const MpsReaderDimensionAndNames& dimension_and_names_;
+  std::vector<Eigen::Triplet<double, int64_t>> triplets_;
+};
+
+}  // namespace
+
+absl::StatusOr<QuadraticProgram> ReadMpsLinearProgram(
+    const std::string& lp_file, bool include_names) {
+  MpsReaderDimensionAndNames dimension_and_names;
+  // Collect MPS format, sizes and names.
+  MPSReaderTemplate<MpsReaderDimensionAndNames> pass_one_reader;
+  const absl::StatusOr<MPSReaderFormat> pass_one_format =
+      pass_one_reader.ParseFile(lp_file, &dimension_and_names,
+                                MPSReaderFormat::kAutoDetect);
+  if (!pass_one_format.ok()) {
+    return util::StatusBuilder(pass_one_format.status()) << absl::StrFormat(
+               "Could not read or parse file `%s` as an MPS file", lp_file);
+  }
+  if (dimension_and_names.FailedToParse()) {
+    return absl::InvalidArgumentError(absl::StrFormat(
+        "Could not read or parse file `%s` as an MPS file, or unsupported "
+        "features/sections found",
+        lp_file));
+  }
+  DCHECK(*pass_one_format == MPSReaderFormat::kFixed ||
+         *pass_one_format == MPSReaderFormat::kFree);
+  // Populate the QP with pre-allocated sizes.
+  MpsReaderQpDataWrapper qp_data_wrapper(&dimension_and_names, include_names);
+  MPSReaderTemplate<MpsReaderQpDataWrapper> pass_two_reader;
+  const absl::StatusOr<MPSReaderFormat> pass_two_format =
+      pass_two_reader.ParseFile(lp_file, &qp_data_wrapper, *pass_one_format);
+  if (!pass_two_format.ok()) {
+    return util::StatusBuilder(pass_two_format.status()) << absl::StrFormat(
+               "Could not read or parse file `%s` as an MPS file "
+               "(maybe file changed between reads?)",
+               lp_file);
+  }
+  DCHECK(*pass_one_format == *pass_two_format);
+  return qp_data_wrapper.GetAndClearQuadraticProgram();
+}
+
+QuadraticProgram ReadMpsLinearProgramOrDie(const std::string& lp_file,
+                                           bool include_names) {
+  const absl::StatusOr<QuadraticProgram> result =
+      ReadMpsLinearProgram(lp_file, include_names);
+  if (!result.ok()) {
+    LOG(QFATAL) << "Error reading MPS Linear Program from " << lp_file << ": "
+                << result.status().message();
+  }
+  return *result;
+}
+
 }  // namespace operations_research::pdlp

ortools/pdlp/quadratic_program_io.h

@@ -20,6 +20,7 @@
 
 #include "absl/container/flat_hash_map.h"
 #include "absl/status/status.h"
+#include "absl/status/statusor.h"
 #include "ortools/pdlp/quadratic_program.h"
 
 namespace operations_research::pdlp {
@@ -34,6 +35,10 @@ QuadraticProgram ReadQuadraticProgramOrDie(const std::string& filename,
 
 QuadraticProgram ReadMpsLinearProgramOrDie(const std::string& lp_file,
                                            bool include_names = false);
+
+absl::StatusOr<QuadraticProgram> ReadMpsLinearProgram(
+    const std::string& lp_file, bool include_names = false);
+
 // The input may be `MPModelProto` in text format, binary format, or JSON,
 // possibly gzipped.
 QuadraticProgram ReadMPModelProtoFileOrDie(