From 3a2510f01c23bbfa5a6731b10cc004ad5d64493e Mon Sep 17 00:00:00 2001 From: Laurent Perron Date: Tue, 11 Jul 2023 15:25:13 -0700 Subject: [PATCH] more pdlp work --- ortools/pdlp/primal_dual_hybrid_gradient.cc | 6 +- ortools/pdlp/quadratic_program_io.cc | 326 +++++++++++++++++++- ortools/pdlp/quadratic_program_io.h | 5 + 3 files changed, 319 insertions(+), 18 deletions(-) diff --git a/ortools/pdlp/primal_dual_hybrid_gradient.cc b/ortools/pdlp/primal_dual_hybrid_gradient.cc index 2c98bbfc9a9..701c081854f 100644 --- a/ortools/pdlp/primal_dual_hybrid_gradient.cc +++ b/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) { diff --git a/ortools/pdlp/quadratic_program_io.cc b/ortools/pdlp/quadratic_program_io.cc index da9414ea82d..5fe4fa7d0c5 100644 --- a/ortools/pdlp/quadratic_program_io.cc +++ b/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 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& ColNameIndexMap() const { + return col_name_to_index_; + } + + // Return a const-hash map of row names to indices. + const absl::flat_hash_map& RowNameIndexMap() const { + return row_name_to_index_; + } + + private: + bool read_or_parse_failed_ = false; + absl::flat_hash_map col_name_to_index_; + absl::flat_hash_map 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 dimension_parser; +// dimension_parser.ParseFile(file_name, &dimension_and_name); +// // Store QP problem coefficients. +// MpsReaderQpDataWrapper qp_wrapper(&dimension_and_name, include_names); +// MPSReaderTemplate 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(dimension_and_names_.NumVariables()); + quadratic_program_.constraint_names = + std::vector(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> triplets_; +}; + +} // namespace + +absl::StatusOr ReadMpsLinearProgram( + const std::string& lp_file, bool include_names) { + MpsReaderDimensionAndNames dimension_and_names; + // Collect MPS format, sizes and names. + MPSReaderTemplate pass_one_reader; + const absl::StatusOr 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 pass_two_reader; + const absl::StatusOr 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 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 diff --git a/ortools/pdlp/quadratic_program_io.h b/ortools/pdlp/quadratic_program_io.h index a122aeeb96a..ea60e9c0cba 100644 --- a/ortools/pdlp/quadratic_program_io.h +++ b/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 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(