From c648ee7b2f21ade8bae2f1f8e3be31b1d42d0d0e Mon Sep 17 00:00:00 2001 From: Peter Mitri Date: Thu, 28 Nov 2024 14:01:41 +0100 Subject: [PATCH 01/14] Add XPRESS solver support to MathOpt --- ortools/linear_solver/xpress_interface.cc | 7 - ortools/math_opt/cpp/parameters.h | 6 + ortools/math_opt/parameters.proto | 6 + ortools/math_opt/solvers/CMakeLists.txt | 28 + ortools/math_opt/solvers/xpress/g_xpress.cc | 235 ++++++++ ortools/math_opt/solvers/xpress/g_xpress.h | 104 ++++ ortools/math_opt/solvers/xpress_solver.cc | 508 ++++++++++++++++++ ortools/math_opt/solvers/xpress_solver.h | 219 ++++++++ .../math_opt/solvers/xpress_solver_test.cc | 85 +++ ortools/xpress/environment.h | 20 + 10 files changed, 1211 insertions(+), 7 deletions(-) create mode 100644 ortools/math_opt/solvers/xpress/g_xpress.cc create mode 100644 ortools/math_opt/solvers/xpress/g_xpress.h create mode 100644 ortools/math_opt/solvers/xpress_solver.cc create mode 100644 ortools/math_opt/solvers/xpress_solver.h create mode 100644 ortools/math_opt/solvers/xpress_solver_test.cc diff --git a/ortools/linear_solver/xpress_interface.cc b/ortools/linear_solver/xpress_interface.cc index 888ba6a4774..9e1abe96697 100644 --- a/ortools/linear_solver/xpress_interface.cc +++ b/ortools/linear_solver/xpress_interface.cc @@ -204,13 +204,6 @@ void interruptXPRESS(XPRSprob& xprsProb, CUSTOM_INTERRUPT_REASON reason) { XPRSinterrupt(xprsProb, 1000 + reason); } -enum XPRS_BASIS_STATUS { - XPRS_AT_LOWER = 0, - XPRS_BASIC = 1, - XPRS_AT_UPPER = 2, - XPRS_FREE_SUPER = 3 -}; - // In case we need to return a double but don't have a value for that // we just return a NaN. #if !defined(XPRS_NAN) diff --git a/ortools/math_opt/cpp/parameters.h b/ortools/math_opt/cpp/parameters.h index ba25d991226..da756cf3ddc 100644 --- a/ortools/math_opt/cpp/parameters.h +++ b/ortools/math_opt/cpp/parameters.h @@ -109,6 +109,12 @@ enum class SolverType { // Slow/not recommended for production. Not an LP solver (no dual information // returned). kSantorini = SOLVER_TYPE_SANTORINI, + + // Fico XPRESS solver (third party). + // + // Supports LP, MIP, and nonconvex integer quadratic problems. + // A fast option, but has special licensing. + kXpress = SOLVER_TYPE_XPRESS }; MATH_OPT_DEFINE_ENUM(SolverType, SOLVER_TYPE_UNSPECIFIED); diff --git a/ortools/math_opt/parameters.proto b/ortools/math_opt/parameters.proto index af0c01bb72c..294595fa28b 100644 --- a/ortools/math_opt/parameters.proto +++ b/ortools/math_opt/parameters.proto @@ -105,6 +105,12 @@ enum SolverTypeProto { // Slow/not recommended for production. Not an LP solver (no dual information // returned). SOLVER_TYPE_SANTORINI = 11; + + // Fico XPRESS solver (third party). + // + // Supports LP, MIP, and nonconvex integer quadratic problems. + // A fast option, but has special licensing. + SOLVER_TYPE_XPRESS = 12; } // Selects an algorithm for solving linear programs. diff --git a/ortools/math_opt/solvers/CMakeLists.txt b/ortools/math_opt/solvers/CMakeLists.txt index e6678ec740f..a10c0ad3cc6 100644 --- a/ortools/math_opt/solvers/CMakeLists.txt +++ b/ortools/math_opt/solvers/CMakeLists.txt @@ -233,3 +233,31 @@ if(USE_HIGHS) "$" ) endif() + +if(USE_XPRESS) + add_subdirectory(xpress) + ortools_cxx_test( + NAME + math_opt_solvers_xpress_solver_test + SOURCES + "xpress_solver_test.cc" + LINK_LIBRARIES + GTest::gmock + GTest::gmock_main + absl::status + ortools::math_opt_matchers + "$" + "$" + "$" + "$" + "$" + #"$" + #"$" + "$" + "$" + "$" + "$" + "$" + "$" + ) +endif() diff --git a/ortools/math_opt/solvers/xpress/g_xpress.cc b/ortools/math_opt/solvers/xpress/g_xpress.cc new file mode 100644 index 00000000000..ee6c8680585 --- /dev/null +++ b/ortools/math_opt/solvers/xpress/g_xpress.cc @@ -0,0 +1,235 @@ +// Copyright 2010-2024 Google LLC +// Licensed under the Apache License, Version 2.0 (the "License"); +// you may not use this file except in compliance with the License. +// You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, software +// distributed under the License is distributed on an "AS IS" BASIS, +// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +// See the License for the specific language governing permissions and +// limitations under the License. + +#include "ortools/math_opt/solvers/xpress/g_xpress.h" + +#include +#include +#include +#include +#include + +#include "absl/log/check.h" +#include "absl/log/die_if_null.h" +#include "absl/memory/memory.h" +#include "absl/status/status.h" +#include "absl/status/statusor.h" +#include "absl/strings/str_format.h" +#include "absl/types/span.h" +#include "ortools/base/logging.h" +#include "ortools/base/source_location.h" +#include "ortools/base/status_builder.h" +#include "ortools/base/status_macros.h" +#include "ortools/xpress/environment.h" + +// The argument to this macro is the invocation of a XPRS function that +// returns a status. If the function returns non-zero the macro aborts +// the program with an appropriate error message. +#define CHECK_STATUS(s) \ + do { \ + int const status_ = s; \ + CHECK_EQ(0, status_); \ + } while (0) + +namespace operations_research::math_opt { +constexpr int kXpressOk = 0; + +absl::Status Xpress::ToStatus(const int xprs_err, + const absl::StatusCode code) const { + if (xprs_err == kXpressOk) { + return absl::OkStatus(); + } + char errmsg[512]; + XPRSgetlasterror(xpress_model_, errmsg); + return util::StatusBuilder(code) + << "Xpress error code: " << xprs_err << ", message: " << errmsg; +} + +Xpress::Xpress(XPRSprob& model) + : xpress_model_(ABSL_DIE_IF_NULL(std::move(model))) {} + +absl::StatusOr> Xpress::New( + const std::string& model_name) { + bool correctlyLoaded = initXpressEnv(); + CHECK(correctlyLoaded); + XPRSprob* model; + CHECK_STATUS(XPRScreateprob(model)); + DCHECK(model != nullptr); // should not be NULL if status=0 + CHECK_STATUS(XPRSaddcbmessage(*model, optimizermsg, NULL, 0)); + return absl::WrapUnique(new Xpress(*model)); +} + +void XPRS_CC optimizermsg(XPRSprob prob, void* data, const char* sMsg, int nLen, + int nMsgLvl) { + printf("%*s\n", nLen, sMsg); +} + +Xpress::~Xpress() { + CHECK_STATUS(XPRSdestroyprob(xpress_model_)); + CHECK_STATUS(XPRSfree()); +} + +absl::Status Xpress::AddVars(const absl::Span obj, + const absl::Span lb, + const absl::Span ub, + const absl::Span vtype) { + return AddVars({}, {}, {}, obj, lb, ub, vtype); +} + +absl::Status Xpress::AddVars(const absl::Span vbegin, + const absl::Span vind, + const absl::Span vval, + const absl::Span obj, + const absl::Span lb, + const absl::Span ub, + const absl::Span vtype) { + CHECK_EQ(vind.size(), vval.size()); + const int num_vars = static_cast(lb.size()); + CHECK_EQ(ub.size(), num_vars); + CHECK_EQ(vtype.size(), num_vars); + double* c_obj = nullptr; + if (!obj.empty()) { + CHECK_EQ(obj.size(), num_vars); + c_obj = const_cast(obj.data()); + } + if (!vbegin.empty()) { + CHECK_EQ(vbegin.size(), num_vars); + } + + CHECK_STATUS(XPRSaddcols(xpress_model_, num_vars, 0, c_obj, nullptr, nullptr, + nullptr, lb.data(), ub.data())); + return absl::OkStatus(); +} + +absl::Status Xpress::AddConstrs(const absl::Span sense, + const absl::Span rhs) { + const int num_cons = static_cast(sense.size()); + CHECK_EQ(rhs.size(), num_cons); + return ToStatus(XPRSaddrows(xpress_model_, num_cons, 0, sense.data(), + rhs.data(), NULL, NULL, NULL, NULL)); +} + +absl::Status Xpress::SetObjective(bool maximize, double offset, + const absl::Span colind, + const absl::Span values) { + int n_cols = colind.size(); + auto c_colind = const_cast(colind.data()); + auto c_values = const_cast(values.data()); + CHECK_STATUS(XPRSchgobj(xpress_model_, n_cols, c_colind, c_values)); + static int indexes[1] = {-1}; + double xprs_values[1] = {-offset}; + CHECK_STATUS(XPRSchgobj(xpress_model_, 1, indexes, xprs_values)); + CHECK_STATUS(XPRSchgobjsense( + xpress_model_, maximize ? XPRS_OBJ_MAXIMIZE : XPRS_OBJ_MINIMIZE)); + return absl::OkStatus(); +} + +absl::Status Xpress::ChgCoeffs(absl::Span rowind, + absl::Span colind, + absl::Span values) { + int n_coefs = rowind.size(); + auto c_rowind = const_cast(rowind.data()); + auto c_colind = const_cast(colind.data()); + auto c_values = const_cast(values.data()); + CHECK_STATUS( + XPRSchgmcoef(xpress_model_, n_coefs, c_rowind, c_colind, c_values)); + return absl::OkStatus(); +} + +absl::StatusOr Xpress::LpOptimizeAndGetStatus() { + RETURN_IF_ERROR(ToStatus(XPRSlpoptimize(xpress_model_, NULL))) + << "XPRESS LP solve failed"; + int xpress_status; + RETURN_IF_ERROR( + ToStatus(XPRSgetintattrib(xpress_model_, XPRS_LPSTATUS, &xpress_status))) + << "Could not get XPRESS status"; + return xpress_status; +} +absl::Status Xpress::PostSolve() { + return ToStatus(XPRSpostsolve(xpress_model_)); +} + +absl::StatusOr Xpress::MipOptimizeAndGetStatus() { + RETURN_IF_ERROR(ToStatus(XPRSmipoptimize(xpress_model_, NULL))) + << "XPRESS MIP solve failed"; + int xpress_status; + RETURN_IF_ERROR( + ToStatus(XPRSgetintattrib(xpress_model_, XPRS_MIPSTATUS, &xpress_status))) + << "Could not get XPRESS status"; + return xpress_status; +} + +void Xpress::Terminate() { XPRSinterrupt(xpress_model_, XPRS_STOP_USER); }; + +absl::StatusOr Xpress::GetIntAttr(int attribute) const { + int result; + RETURN_IF_ERROR(ToStatus(XPRSgetintattrib(xpress_model_, attribute, &result))) + << "Error getting Xpress int attribute: " << attribute; + return result; +} + +absl::Status Xpress::SetIntAttr(int attribute, int value) { + return ToStatus(XPRSsetintcontrol(xpress_model_, attribute, value)); +} + +absl::StatusOr Xpress::GetDoubleAttr(int attribute) const { + double result; + RETURN_IF_ERROR(ToStatus(XPRSgetdblattrib(xpress_model_, attribute, &result))) + << "Error getting Xpress double attribute: " << attribute; + return result; +} + +absl::StatusOr> Xpress::GetPrimalValues() const { + int nCols = GetNumberOfColumns(); + XPRSgetintattrib(xpress_model_, XPRS_COLS, &nCols); + double values[nCols]; + RETURN_IF_ERROR( + ToStatus(XPRSgetlpsol(xpress_model_, values, nullptr, nullptr, nullptr))) + << "Error getting Xpress LP solution"; + std::vector result(values, values + nCols); + return std::move(result); +} + +int Xpress::GetNumberOfRows() const { + int n; + XPRSgetintattrib(xpress_model_, XPRS_ROWS, &n); + return n; +} + +int Xpress::GetNumberOfColumns() const { + int n; + XPRSgetintattrib(xpress_model_, XPRS_COLS, &n); + return n; +} + +std::vector Xpress::GetConstraintDuals() const { + int nRows = GetNumberOfRows(); + double values[nRows]; + CHECK_STATUS(XPRSgetlpsol(xpress_model_, nullptr, nullptr, values, nullptr)); + return {values, values + nRows}; +} +std::vector Xpress::GetReducedCostValues() const { + int nCols = GetNumberOfColumns(); + double values[nCols]; + CHECK_STATUS(XPRSgetlpsol(xpress_model_, nullptr, nullptr, nullptr, values)); + return {values, values + nCols}; +} + +std::vector Xpress::GetVariableBasis() const { + int const nCols = GetNumberOfColumns(); + std::vector basis(nCols); + CHECK_STATUS(XPRSgetbasis(xpress_model_, basis.data(), 0)); + return basis; +} + +} // namespace operations_research::math_opt \ No newline at end of file diff --git a/ortools/math_opt/solvers/xpress/g_xpress.h b/ortools/math_opt/solvers/xpress/g_xpress.h new file mode 100644 index 00000000000..9d5fee2909b --- /dev/null +++ b/ortools/math_opt/solvers/xpress/g_xpress.h @@ -0,0 +1,104 @@ +// Copyright 2010-2024 Google LLC +// Licensed under the Apache License, Version 2.0 (the "License"); +// you may not use this file except in compliance with the License. +// You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, software +// distributed under the License is distributed on an "AS IS" BASIS, +// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +// See the License for the specific language governing permissions and +// limitations under the License. + +// Google C++ bindings for Xpress C API. +// +// Attempts to be as close to the Xpress C API as possible, with the following +// differences: +// * Use destructors to automatically clean up the environment and model. +// * Use absl::Status to propagate errors. +// * Use absl::StatusOr instead of output arguments. +// * Use absl::Span instead of T* and size for array args. +// * Use std::string instead of null terminated char* for string values (note +// that attribute names are still char*). +// * When setting array data, accept const data (absl::Span). +#ifndef OR_TOOLS_MATH_OPT_SOLVERS_XPRESS_G_XPRESS_H_ +#define OR_TOOLS_MATH_OPT_SOLVERS_XPRESS_G_XPRESS_H_ + +#include +#include +#include +#include +#include + +#include "absl/status/status.h" +#include "absl/status/statusor.h" +#include "absl/types/span.h" +#include "ortools/xpress/environment.h" + +namespace operations_research::math_opt { + +class Xpress { + public: + Xpress() = delete; + + absl::Status ToStatus( + int xprs_err, + absl::StatusCode code = absl::StatusCode::kInvalidArgument) const; + + // Creates a new Xpress + static absl::StatusOr> New( + const std::string& model_name); + + ~Xpress(); + + absl::StatusOr GetIntAttr(int attribute) const; + absl::Status SetIntAttr(int attribute, int value); + + absl::StatusOr GetDoubleAttr(int attribute) const; + + absl::Status AddVars(absl::Span obj, + absl::Span lb, absl::Span ub, + absl::Span vtype); + + absl::Status AddVars(absl::Span vbegin, absl::Span vind, + absl::Span vval, + absl::Span obj, + absl::Span lb, absl::Span ub, + absl::Span vtype); + + absl::Status AddConstrs(absl::Span sense, + absl::Span rhs); + absl::Status SetObjective(bool maximize, double offset, + absl::Span colind, + absl::Span values); + + absl::Status ChgCoeffs(absl::Span cind, absl::Span vind, + absl::Span val); + + absl::StatusOr LpOptimizeAndGetStatus(); + absl::StatusOr MipOptimizeAndGetStatus(); + absl::Status PostSolve(); + + void Terminate(); + + absl::StatusOr> GetPrimalValues() const; + std::vector GetConstraintDuals() const; + std::vector GetReducedCostValues() const; + std::vector GetVariableBasis() const; + + private: + XPRSprob xpress_model_; + + explicit Xpress(XPRSprob& model); + + int GetNumberOfRows() const; + int GetNumberOfColumns() const; + + static void XPRS_CC optimizermsg(XPRSprob prob, void* data, const char* sMsg, + int nLen, int nMsgLvl); +}; + +} // namespace operations_research::math_opt + +#endif // OR_TOOLS_MATH_OPT_SOLVERS_XPRESS_G_XPRESS_H_ diff --git a/ortools/math_opt/solvers/xpress_solver.cc b/ortools/math_opt/solvers/xpress_solver.cc new file mode 100644 index 00000000000..71e92d72d93 --- /dev/null +++ b/ortools/math_opt/solvers/xpress_solver.cc @@ -0,0 +1,508 @@ +// Copyright 2010-2024 Google LLC +// Licensed under the Apache License, Version 2.0 (the "License"); +// you may not use this file except in compliance with the License. +// You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, software +// distributed under the License is distributed on an "AS IS" BASIS, +// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +// See the License for the specific language governing permissions and +// limitations under the License. + +#include "xpress_solver.h" + +#include "ortools/base/map_util.h" +#include "ortools/base/protoutil.h" +#include "ortools/base/status_macros.h" +#include "ortools/math_opt/core/math_opt_proto_utils.h" +#include "ortools/math_opt/cpp/solve_result.h" +#include "ortools/xpress/environment.h" + +namespace operations_research::math_opt { + +constexpr SupportedProblemStructures kXpressSupportedStructures = { + .integer_variables = SupportType::kNotImplemented, + .multi_objectives = SupportType::kNotImplemented, + .quadratic_objectives = SupportType::kNotImplemented, + .quadratic_constraints = SupportType::kNotImplemented, + .second_order_cone_constraints = SupportType::kNotImplemented, + .sos1_constraints = SupportType::kNotImplemented, + .sos2_constraints = SupportType::kNotImplemented, + .indicator_constraints = SupportType::kNotImplemented}; + +absl::StatusOr> XpressSolver::New( + const ModelProto& input_model, const SolverInterface::InitArgs& init_args) { + if (!XpressIsCorrectlyInstalled()) { + return absl::InvalidArgumentError("Xpress is not correctly installed."); + } + RETURN_IF_ERROR( + ModelIsSupported(input_model, kXpressSupportedStructures, "XPRESS")); + + // We can add here extra checks that are not made in ModelIsSupported + // (for example, if XPRESS does not support multi-objective with quad terms) + + ASSIGN_OR_RETURN(std::unique_ptr xpr, + Xpress::New(input_model.name())); + auto xpress_solver = absl::WrapUnique(new XpressSolver(std::move(xpr))); + RETURN_IF_ERROR(xpress_solver->LoadModel(input_model)); + return xpress_solver; +} + +absl::Status XpressSolver::LoadModel(const ModelProto& input_model) { + CHECK(xpress_ != nullptr); + // TODO: set prob name, use XPRSsetprobname (must be added to environment) + // (!) must be truncated to MAXPROBNAMELENGTH + // RETURN_IF_ERROR(xpress_->SetProbName(input_model.name()); + RETURN_IF_ERROR(AddNewVariables(input_model.variables())); + RETURN_IF_ERROR(AddNewLinearConstraints(input_model.linear_constraints())); + // TODO: instead of changing coefficients, set them when adding constraints? + RETURN_IF_ERROR(ChangeCoefficients(input_model.linear_constraint_matrix())); + RETURN_IF_ERROR(AddSingleObjective(input_model.objective())); + return absl::OkStatus(); +} +absl::Status XpressSolver::AddNewVariables( + const VariablesProto& new_variables) { + const int num_new_variables = new_variables.lower_bounds().size(); + std::vector variable_type(num_new_variables); + for (int j = 0; j < num_new_variables; ++j) { + const VariableId id = new_variables.ids(j); + InsertOrDie(&variables_map_, id, j + num_xpress_variables_); + variable_type[j] = + new_variables.integers(j) ? XPRS_INTEGER : XPRS_CONTINUOUS; + if (new_variables.integers(j)) { + is_mip_ = true; + return absl::UnimplementedError("XpressSolver does not handle MIPs yet"); + } + } + RETURN_IF_ERROR(xpress_->AddVars({}, new_variables.lower_bounds(), + new_variables.upper_bounds(), + variable_type)); + num_xpress_variables_ += num_new_variables; + + // Not adding names for performance (have to call XPRSaddnames) + // TODO : keep names in a cache and add them when needed? + + return absl::OkStatus(); +} + +XpressSolver::XpressSolver(std::unique_ptr g_xpress) + : xpress_(std::move(g_xpress)) {} + +absl::Status XpressSolver::AddNewLinearConstraints( + const LinearConstraintsProto& constraints) { + const int num_new_constraints = constraints.lower_bounds().size(); + + // Constraints are translated into: + // 1. ax <= upper_bound (if lower bound <= -INFINITY, and upper_bound + // is finite and less than INFINITY) + // 2. ax >= lower_bound (if upper bound >= INFINITY, and lower_bound is + // finite and greater than -INFINITY) + // 3. ax == xxxxx_bound (if both bounds are finite, equal, and their + // absolute values less than INFINITY) + // 4. ax - slack = 0.0 (otherwise, + // slack bounds == [lower_bound, upper_bound]) + std::vector constraint_rhs; + std::vector constraint_sense; + std::vector new_slacks; + constraint_rhs.reserve(num_new_constraints); + constraint_sense.reserve(num_new_constraints); + new_slacks.reserve(num_new_constraints); + for (int i = 0; i < num_new_constraints; ++i) { + const int64_t id = constraints.ids(i); + LinearConstraintData& constraint_data = + gtl::InsertKeyOrDie(&linear_constraints_map_, id); + const double lb = constraints.lower_bounds(i); + const double ub = constraints.upper_bounds(i); + constraint_data.lower_bound = lb; + constraint_data.upper_bound = ub; + constraint_data.constraint_index = i + num_xpress_lin_cons_; + char sense = XPRS_EQUAL; + double rhs = 0.0; + const bool lb_is_xprs_neg_inf = lb <= XPRS_MINUSINFINITY; + const bool ub_is_xprs_pos_inf = ub >= XPRS_PLUSINFINITY; + if (lb_is_xprs_neg_inf && !ub_is_xprs_pos_inf) { + sense = XPRS_LESS_EQUAL; + rhs = ub; + } else if (!lb_is_xprs_neg_inf && ub_is_xprs_pos_inf) { + sense = XPRS_GREATER_EQUAL; + rhs = lb; + } else if (lb == ub) { + sense = XPRS_EQUAL; + rhs = lb; + } else { + // Note that constraints where the lower bound and the upper bound are + // -+infinity translate into a range constraint with an unbounded slack. + constraint_data.slack_index = new_slacks.size() + num_xpress_variables_; + new_slacks.push_back(&constraint_data); + } + constraint_rhs.emplace_back(rhs); + constraint_sense.emplace_back(sense); + } + // Add all constraints in one call. + RETURN_IF_ERROR(xpress_->AddConstrs(constraint_sense, constraint_rhs)); + num_xpress_lin_cons_ += num_new_constraints; + // Add slacks for true ranged constraints (if needed) + if (!new_slacks.empty()) { + // TODO + return absl::UnimplementedError( + "XpressSolver does not support slack variables yet"); + } + return absl::OkStatus(); +} + +absl::Status XpressSolver::AddSingleObjective(const ObjectiveProto& objective) { + /*if (objective.has_quadratic_coefficients()) { + return absl::UnimplementedError( + "Quadratic objectives are not yet implemented in XPRESS solver " + "interface."); + }*/ + if (objective.linear_coefficients().ids_size() == 0) { + return absl::OkStatus(); + } + std::vector index; + index.reserve(objective.linear_coefficients().ids_size()); + for (const int64_t id : objective.linear_coefficients().ids()) { + index.push_back(variables_map_.at(id)); + } + RETURN_IF_ERROR( + xpress_->SetObjective(objective.maximize(), objective.offset(), index, + objective.linear_coefficients().values())); + is_maximize_ = objective.maximize(); + return absl::OkStatus(); +} + +absl::Status XpressSolver::ChangeCoefficients( + const SparseDoubleMatrixProto& matrix) { + const int num_coefficients = matrix.row_ids().size(); + std::vector row_index(num_coefficients); + std::vector col_index(num_coefficients); + for (int k = 0; k < num_coefficients; ++k) { + row_index[k] = + linear_constraints_map_.at(matrix.row_ids(k)).constraint_index; + col_index[k] = variables_map_.at(matrix.column_ids(k)); + } + return xpress_->ChgCoeffs(row_index, col_index, matrix.coefficients()); +} + +absl::StatusOr XpressSolver::Solve( + const SolveParametersProto& parameters, + const ModelSolveParametersProto& model_parameters, + MessageCallback message_cb, + const CallbackRegistrationProto& callback_registration, Callback cb, + const SolveInterrupter* interrupter) { + const absl::Time start = absl::Now(); + + // TODO: set solve parameters + // TODO: set basis + // TODO: set hints + // TODO: set branching properties + // TODO: set lazy constraints + // TODO: add interrupter using xpress_->Terminate(); + + RETURN_IF_ERROR(CallXpressSolve(parameters.enable_output())) + << "Error during XPRESS solve"; + + ASSIGN_OR_RETURN(SolveResultProto solve_result, + ExtractSolveResultProto(start, model_parameters)); + + return solve_result; +} + +absl::Status XpressSolver::CallXpressSolve(bool enableOutput) { + // Enable screen output right before solve + if (enableOutput) { + RETURN_IF_ERROR(xpress_->SetIntAttr(XPRS_OUTPUTLOG, 1)) + << "Unable to enable XPRESS logs"; + } + // Solve + if (is_mip_) { + ASSIGN_OR_RETURN(xpress_status_, xpress_->MipOptimizeAndGetStatus()); + } else { + ASSIGN_OR_RETURN(xpress_status_, xpress_->LpOptimizeAndGetStatus()); + } + // Post-solve + if (!(is_mip_ ? (xpress_status_ == XPRS_MIP_OPTIMAL) + : (xpress_status_ == XPRS_LP_OPTIMAL))) { + RETURN_IF_ERROR(xpress_->PostSolve()) << "Post-solve failed in XPRESS"; + } + // Disable screen output right after solve + if (enableOutput) { + RETURN_IF_ERROR(xpress_->SetIntAttr(XPRS_OUTPUTLOG, 0)) + << "Unable to disable XPRESS logs"; + } + return absl::OkStatus(); +} + +absl::StatusOr XpressSolver::ExtractSolveResultProto( + absl::Time start, const ModelSolveParametersProto& model_parameters) { + SolveResultProto result; + ASSIGN_OR_RETURN(SolutionsAndClaims solution_and_claims, + GetSolutions(model_parameters)); + for (SolutionProto& solution : solution_and_claims.solutions) { + *result.add_solutions() = std::move(solution); + } + ASSIGN_OR_RETURN(*result.mutable_solve_stats(), GetSolveStats(start)); + ASSIGN_OR_RETURN(const double best_primal_bound, GetBestPrimalBound()); + ASSIGN_OR_RETURN(const double best_dual_bound, GetBestDualBound()); + auto solution_claims = solution_and_claims.solution_claims; + ASSIGN_OR_RETURN(*result.mutable_termination(), + ConvertTerminationReason(solution_claims, best_primal_bound, + best_dual_bound)); + return std::move(result); +} + +absl::StatusOr XpressSolver::GetBestPrimalBound() const { + ASSIGN_OR_RETURN( + const double best_primal_bound, + xpress_->GetDoubleAttr(is_mip_ ? XPRS_MIPOBJVAL : XPRS_LPOBJVAL)); + return best_primal_bound; +} + +absl::StatusOr XpressSolver::GetBestDualBound() const { + ASSIGN_OR_RETURN(const double best_dual_bound, + xpress_->GetDoubleAttr(XPRS_BESTBOUND)); + return best_dual_bound; +} + +absl::StatusOr XpressSolver::GetSolutions( + const ModelSolveParametersProto& model_parameters) { + if (is_mip_) { + return absl::UnimplementedError("XpressSolver does not handle MIPs yet"); + } else { + return GetLpSolution(model_parameters); + } +} + +absl::StatusOr XpressSolver::GetLpSolution( + const ModelSolveParametersProto& model_parameters) { + ASSIGN_OR_RETURN(auto primal_solution_and_claim, + GetConvexPrimalSolutionIfAvailable(model_parameters)); + ASSIGN_OR_RETURN(auto dual_solution_and_claim, + GetConvexDualSolutionIfAvailable(model_parameters)); + ASSIGN_OR_RETURN(auto basis, GetBasisIfAvailable()); + const SolutionClaims solution_claims = { + .primal_feasible_solution_exists = + primal_solution_and_claim.feasible_solution_exists, + .dual_feasible_solution_exists = + dual_solution_and_claim.feasible_solution_exists}; + + if (!primal_solution_and_claim.solution.has_value() && + !dual_solution_and_claim.solution.has_value() && !basis.has_value()) { + return SolutionsAndClaims{.solution_claims = solution_claims}; + } + SolutionsAndClaims solution_and_claims{.solution_claims = solution_claims}; + SolutionProto& solution = + solution_and_claims.solutions.emplace_back(SolutionProto()); + if (primal_solution_and_claim.solution.has_value()) { + *solution.mutable_primal_solution() = + std::move(*primal_solution_and_claim.solution); + } + if (dual_solution_and_claim.solution.has_value()) { + *solution.mutable_dual_solution() = + std::move(*dual_solution_and_claim.solution); + } + if (basis.has_value()) { + *solution.mutable_basis() = std::move(*basis); + } + return solution_and_claims; +} + +bool XpressSolver::isFeasible() const { + return xpress_status_ == (is_mip_ ? XPRS_MIP_OPTIMAL : XPRS_LP_OPTIMAL); +} + +SolutionStatusProto XpressSolver::getLpSolutionStatus() const { + switch (xpress_status_) { + case XPRS_LP_OPTIMAL: + return SOLUTION_STATUS_FEASIBLE; + case XPRS_LP_INFEAS: + return SOLUTION_STATUS_INFEASIBLE; + case XPRS_LP_UNBOUNDED: + return SOLUTION_STATUS_UNDETERMINED; + default: + return SOLUTION_STATUS_UNSPECIFIED; + } +} + +absl::StatusOr> +XpressSolver::GetConvexPrimalSolutionIfAvailable( + const ModelSolveParametersProto& model_parameters) const { + PrimalSolutionProto primal_solution; + if (isFeasible()) { + ASSIGN_OR_RETURN(const double sol_val, + xpress_->GetDoubleAttr(XPRS_LPOBJVAL)); + primal_solution.set_objective_value(sol_val); + } else { + // TODO + } + primal_solution.set_feasibility_status(getLpSolutionStatus()); + XpressVectorToSparseDoubleVector(xpress_->GetPrimalValues().value(), + variables_map_, + *primal_solution.mutable_variable_values(), + model_parameters.variable_values_filter()); + const bool primal_feasible_solution_exists = + (primal_solution.feasibility_status() == SOLUTION_STATUS_FEASIBLE); + return SolutionAndClaim{ + .solution = std::move(primal_solution), + .feasible_solution_exists = primal_feasible_solution_exists}; +} + +absl::StatusOr> +XpressSolver::GetConvexDualSolutionIfAvailable( + const ModelSolveParametersProto& model_parameters) const { + DualSolutionProto dual_solution; + const std::vector xprs_constraint_duals = + xpress_->GetConstraintDuals(); + XpressVectorToSparseDoubleVector(xprs_constraint_duals, + linear_constraints_map_, + *dual_solution.mutable_dual_values(), + model_parameters.dual_values_filter()); + + const std::vector xprs_reduced_cost_values = + xpress_->GetReducedCostValues(); + XpressVectorToSparseDoubleVector(xprs_reduced_cost_values, variables_map_, + *dual_solution.mutable_reduced_costs(), + model_parameters.reduced_costs_filter()); + + if (isFeasible()) { + ASSIGN_OR_RETURN(const double sol_val, + xpress_->GetDoubleAttr(XPRS_LPOBJVAL)); + dual_solution.set_objective_value(sol_val); + } else { + // TODO + } + dual_solution.set_feasibility_status(getLpSolutionStatus()); + bool dual_feasible_solution_exists = + (dual_solution.feasibility_status() == SOLUTION_STATUS_FEASIBLE); + ASSIGN_OR_RETURN(const double best_dual_bound, GetBestDualBound()); + if (dual_feasible_solution_exists || std::isfinite(best_dual_bound)) { + dual_feasible_solution_exists = true; + } else if (xpress_status_ == XPRS_LP_OPTIMAL) { + return absl::InternalError( + "Xpress status is XPRS_LP_OPTIMAL, but XPRS_BESTBOUND is " + "unavailable or infinite, and no dual feasible solution is returned"); + } + return SolutionAndClaim{ + .solution = std::move(dual_solution), + .feasible_solution_exists = dual_feasible_solution_exists}; +} + +inline BasisStatusProto ConvertVariableStatus(const int status) { + switch (status) { + case XPRS_BASIC: + return BASIS_STATUS_BASIC; + case XPRS_AT_LOWER: + return BASIS_STATUS_AT_LOWER_BOUND; + case XPRS_AT_UPPER: + return BASIS_STATUS_AT_UPPER_BOUND; + case XPRS_FREE_SUPER: + return BASIS_STATUS_FREE; + default: + return BASIS_STATUS_UNSPECIFIED; + } +} + +absl::StatusOr> XpressSolver::GetBasisIfAvailable() { + // Variable basis + BasisProto basis; + auto xprs_variable_basis_status = xpress_->GetVariableBasis(); + for (auto [variable_id, xprs_variable_index] : variables_map_) { + basis.mutable_variable_status()->add_ids(variable_id); + const BasisStatusProto variable_status = + ConvertVariableStatus(xprs_variable_basis_status[xprs_variable_index]); + if (variable_status == BASIS_STATUS_UNSPECIFIED) { + return absl::InternalError( + absl::StrCat("Invalid Xpress variable basis status: ", + xprs_variable_basis_status[xprs_variable_index])); + } + basis.mutable_variable_status()->add_values(variable_status); + } + // Constraint basis + // TODO : implement this, mocked for now (else Basis validation will fail) + for (auto [constraint_id, xprs_ct_index] : linear_constraints_map_) { + basis.mutable_constraint_status()->add_ids(constraint_id); + basis.mutable_constraint_status()->add_values(BASIS_STATUS_BASIC); + } + // Dual basis + basis.set_basic_dual_feasibility(SOLUTION_STATUS_UNDETERMINED); + if (xpress_status_ == XPRS_LP_OPTIMAL) { + basis.set_basic_dual_feasibility(SOLUTION_STATUS_FEASIBLE); + } else if (xpress_status_ == XPRS_LP_UNBOUNDED) { + basis.set_basic_dual_feasibility(SOLUTION_STATUS_INFEASIBLE); + } + return std::move(basis); +} + +absl::StatusOr XpressSolver::GetSolveStats( + absl::Time start) const { + SolveStatsProto solve_stats; + CHECK_OK(util_time::EncodeGoogleApiProto(absl::Now() - start, + solve_stats.mutable_solve_time())); + // TODO : complete these stats + return solve_stats; +} + +template +void XpressSolver::XpressVectorToSparseDoubleVector( + absl::Span xpress_values, const T& map, + SparseDoubleVectorProto& result, + const SparseVectorFilterProto& filter) const { + SparseVectorFilterPredicate predicate(filter); + for (auto [id, xpress_data] : map) { + const double value = xpress_values[get_model_index(xpress_data)]; + if (predicate.AcceptsAndUpdate(id, value)) { + result.add_ids(id); + result.add_values(value); + } + } +} + +absl::StatusOr XpressSolver::ConvertTerminationReason( + SolutionClaims solution_claims, double best_primal_bound, + double best_dual_bound) const { + // TODO : improve this + if (!is_mip_) { + switch (xpress_status_) { + case XPRS_LP_OPTIMAL: + return OptimalTerminationProto(best_primal_bound, best_dual_bound); + case XPRS_LP_INFEAS: + return InfeasibleTerminationProto( + is_maximize_, solution_claims.dual_feasible_solution_exists + ? FEASIBILITY_STATUS_FEASIBLE + : FEASIBILITY_STATUS_UNDETERMINED); + case XPRS_LP_UNBOUNDED: + if (solution_claims.primal_feasible_solution_exists) { + return UnboundedTerminationProto(is_maximize_); + } + return InfeasibleOrUnboundedTerminationProto( + is_maximize_, FEASIBILITY_STATUS_INFEASIBLE, + "Xpress status XPRS_LP_UNBOUNDED"); + default: + return TerminateForReason(is_maximize_, TERMINATION_REASON_IMPRECISE); + } + } else { + return absl::UnimplementedError("XpressSolver does not handle MIPs yet"); + } +} + +absl::StatusOr XpressSolver::Update( + const ModelUpdateProto& model_update) { + // TODO: implement this + return absl::UnimplementedError( + "XpressSolver::Update is not implemented yet"); +} + +absl::StatusOr +XpressSolver::ComputeInfeasibleSubsystem(const SolveParametersProto& parameters, + MessageCallback message_cb, + const SolveInterrupter* interrupter) { + // TODO: implement this + return absl::UnimplementedError( + "XpressSolver::ComputeInfeasibleSubsystem is not implemented yet"); +} + +MATH_OPT_REGISTER_SOLVER(SOLVER_TYPE_XPRESS, XpressSolver::New) + +} // namespace operations_research::math_opt diff --git a/ortools/math_opt/solvers/xpress_solver.h b/ortools/math_opt/solvers/xpress_solver.h new file mode 100644 index 00000000000..4b250c09f51 --- /dev/null +++ b/ortools/math_opt/solvers/xpress_solver.h @@ -0,0 +1,219 @@ +// Copyright 2010-2024 Google LLC +// Licensed under the Apache License, Version 2.0 (the "License"); +// you may not use this file except in compliance with the License. +// You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, software +// distributed under the License is distributed on an "AS IS" BASIS, +// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +// See the License for the specific language governing permissions and +// limitations under the License. + +#ifndef OR_TOOLS_MATH_OPT_SOLVERS_XPRESS_SOLVER_H_ +#define OR_TOOLS_MATH_OPT_SOLVERS_XPRESS_SOLVER_H_ + +#include +#include +#include +#include +#include + +#include "absl/status/status.h" +#include "absl/status/statusor.h" +#include "absl/time/time.h" +#include "absl/types/span.h" +#include "ortools/base/linked_hash_map.h" +#include "ortools/math_opt/callback.pb.h" +#include "ortools/math_opt/core/solver_interface.h" +#include "ortools/math_opt/infeasible_subsystem.pb.h" +#include "ortools/math_opt/model.pb.h" +#include "ortools/math_opt/model_parameters.pb.h" +#include "ortools/math_opt/model_update.pb.h" +#include "ortools/math_opt/parameters.pb.h" +#include "ortools/math_opt/result.pb.h" +#include "ortools/math_opt/solution.pb.h" +#include "ortools/math_opt/solvers/xpress/g_xpress.h" +#include "ortools/math_opt/sparse_containers.pb.h" +#include "ortools/util/solve_interrupter.h" + +namespace operations_research::math_opt { + +// Interface to FICO XPRESS solver +// Largely inspired by the Gurobi interface +class XpressSolver : public SolverInterface { + public: + // Creates the XPRESS solver and loads the model into it + static absl::StatusOr> New( + const ModelProto& input_model, + const SolverInterface::InitArgs& init_args); + + // Solves the optimization problem + absl::StatusOr Solve( + const SolveParametersProto& parameters, + const ModelSolveParametersProto& model_parameters, + MessageCallback message_cb, + const CallbackRegistrationProto& callback_registration, Callback cb, + const SolveInterrupter* interrupter) override; + absl::Status CallXpressSolve(bool enableOutput); + + // Updates the problem (not implemented yet) + absl::StatusOr Update(const ModelUpdateProto& model_update) override; + + // Computes the infeasible subsystem (not implemented yet) + absl::StatusOr + ComputeInfeasibleSubsystem(const SolveParametersProto& parameters, + MessageCallback message_cb, + const SolveInterrupter* interrupter) override; + + private: + explicit XpressSolver(std::unique_ptr g_xpress); + + // For easing reading the code, we declare these types: + using VariableId = int64_t; + using AuxiliaryObjectiveId = int64_t; + using LinearConstraintId = int64_t; + using QuadraticConstraintId = int64_t; + using SecondOrderConeConstraintId = int64_t; + using Sos1ConstraintId = int64_t; + using Sos2ConstraintId = int64_t; + using IndicatorConstraintId = int64_t; + using AnyConstraintId = int64_t; + using XpressVariableIndex = int; + using XpressMultiObjectiveIndex = int; + using XpressLinearConstraintIndex = int; + using XpressQuadraticConstraintIndex = int; + using XpressSosConstraintIndex = int; + // A collection of other constraints (e.g., norm, max, indicator) supported by + // Xpress. All general constraints share the same index set. See for more + // detail: https://www.xpress.com/documentation/9.5/refman/constraints.html. + using XpressGeneralConstraintIndex = int; + using XpressAnyConstraintIndex = int; + + static constexpr XpressVariableIndex kUnspecifiedIndex = -1; + static constexpr XpressAnyConstraintIndex kUnspecifiedConstraint = -2; + static constexpr double kPlusInf = XPRS_PLUSINFINITY; + static constexpr double kMinusInf = XPRS_MINUSINFINITY; + + // Data associated with each linear constraint. With it we know if the + // underlying representation is either: + // linear_terms <= upper_bound (if lower bound <= -INFINITY) + // linear_terms >= lower_bound (if upper bound >= INFINTY) + // linear_terms == xxxxx_bound (if upper_bound == lower_bound) + // linear_term - slack == 0 (with slack bounds equal to xxxxx_bound) + struct LinearConstraintData { + XpressLinearConstraintIndex constraint_index = kUnspecifiedConstraint; + // only valid for true ranged constraints. + XpressVariableIndex slack_index = kUnspecifiedIndex; + double lower_bound = kMinusInf; + double upper_bound = kPlusInf; + }; + + struct SolutionClaims { + bool primal_feasible_solution_exists; + bool dual_feasible_solution_exists; + }; + + struct SolutionsAndClaims { + std::vector solutions; + SolutionClaims solution_claims; + }; + + template + struct SolutionAndClaim { + std::optional solution; + bool feasible_solution_exists = false; + }; + + absl::StatusOr ExtractSolveResultProto( + absl::Time start, const ModelSolveParametersProto& model_parameters); + absl::StatusOr GetSolutions( + const ModelSolveParametersProto& model_parameters); + absl::StatusOr GetSolveStats(absl::Time start) const; + + absl::StatusOr GetBestDualBound() const; + absl::StatusOr GetBestPrimalBound() const; + + absl::StatusOr ConvertTerminationReason( + SolutionClaims solution_claims, double best_primal_bound, + double best_dual_bound) const; + + // Returns solution information appropriate and available for an LP (linear + // constraints + linear objective, only). + absl::StatusOr GetLpSolution( + const ModelSolveParametersProto& model_parameters); + bool isFeasible() const; + + // return bool field should be true if a primal solution exists. + absl::StatusOr> + GetConvexPrimalSolutionIfAvailable( + const ModelSolveParametersProto& model_parameters) const; + absl::StatusOr> + GetConvexDualSolutionIfAvailable( + const ModelSolveParametersProto& model_parameters) const; + absl::StatusOr> GetBasisIfAvailable(); + + absl::Status AddNewLinearConstraints( + const LinearConstraintsProto& constraints); + absl::Status AddNewVariables(const VariablesProto& new_variables); + absl::Status AddSingleObjective(const ObjectiveProto& objective); + absl::Status ChangeCoefficients(const SparseDoubleMatrixProto& matrix); + + absl::Status LoadModel(const ModelProto& input_model); + + // Fills in result with the values in xpress_values aided by the index + // conversion from map which should be either variables_map_ or + // linear_constraints_map_ as appropriate. Only key/value pairs that passes + // the filter predicate are added. + template + void XpressVectorToSparseDoubleVector( + absl::Span xpress_values, const T& map, + SparseDoubleVectorProto& result, + const SparseVectorFilterProto& filter) const; + + const std::unique_ptr xpress_; + + // Note that we use linked_hash_map for the indices of the xpress_model_ + // variables and linear constraints to ensure that iteration over the map + // maintains their insertion order (and, thus, the order in which they appear + // in the model). + + // Internal correspondence from variable proto IDs to Xpress-numbered + // variables. + gtl::linked_hash_map variables_map_; + // Internal correspondence from linear constraint proto IDs to + // Xpress-numbered linear constraint and extra information. + gtl::linked_hash_map + linear_constraints_map_; + + int get_model_index(XpressVariableIndex index) const { return index; } + int get_model_index(const LinearConstraintData& index) const { + return index.constraint_index; + } + + SolutionStatusProto getLpSolutionStatus() const; + + // Fields to track the number of Xpress variables and constraints. These + // quantities are updated immediately after adding or removing to the model, + // so it is correct even if XPRESS C API has not yet been called. + + // Number of Xpress variables. + int num_xpress_variables_ = 0; + // Number of Xpress linear constraints. + int num_xpress_lin_cons_ = 0; + // Number of Xpress quadratic constraints. + int num_xpress_quad_cons_ = 0; + // Number of Xpress SOS constraints. + int num_xpress_sos_cons_ = 0; + // Number of Xpress general constraints. + int num_xpress_gen_cons_ = 0; + + bool is_mip_ = false; + bool is_maximize_ = false; + int xpress_status_ = 0; +}; + +} // namespace operations_research::math_opt + +#endif // OR_TOOLS_MATH_OPT_SOLVERS_XPRESS_SOLVER_H_ diff --git a/ortools/math_opt/solvers/xpress_solver_test.cc b/ortools/math_opt/solvers/xpress_solver_test.cc new file mode 100644 index 00000000000..87eadab6ab5 --- /dev/null +++ b/ortools/math_opt/solvers/xpress_solver_test.cc @@ -0,0 +1,85 @@ +// Copyright 2010-2024 Google LLC +// Licensed under the Apache License, Version 2.0 (the "License"); +// you may not use this file except in compliance with the License. +// You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, software +// distributed under the License is distributed on an "AS IS" BASIS, +// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +// See the License for the specific language governing permissions and +// limitations under the License. + +#include + +#include + +#include "absl/log/check.h" +#include "absl/status/statusor.h" +#include "gtest/gtest.h" +#include "ortools/base/gmock.h" +#include "ortools/base/init_google.h" +#include "ortools/math_opt/core/solver.h" +#include "ortools/math_opt/cpp/matchers.h" +#include "ortools/math_opt/cpp/math_opt.h" +#include "ortools/math_opt/solver_tests/lp_model_solve_parameters_tests.h" +#include "ortools/math_opt/solver_tests/lp_tests.h" +#include "ortools/port/proto_utils.h" + +using namespace operations_research::math_opt; + +using ::testing::AnyOf; +using ::testing::status::IsOkAndHolds; +constexpr double kTolerance = 1.0e-5; +constexpr double kInf = std::numeric_limits::infinity(); + +struct SolvedModel { + std::unique_ptr model; + SolveResult expected_result; +}; + +// max(x + 2 * y) +// s.t.: -1 <= x <= 1.5 +// 0 <= y <= 1 +// x + y <= 1.5 +// optimal solution is: (x, y) = (0.5, 1) +// obj_value = 2.5 +TEST(XpressSolverTest, SimpleLpTwoVar) { + // Build the model. + Model lp_model("getting_started_lp"); + const Variable x = lp_model.AddContinuousVariable(-1.0, 1.5, "x"); + const Variable y = lp_model.AddContinuousVariable(0.0, 1.0, "y"); + lp_model.AddLinearConstraint(x + y <= 1.5, "c"); + lp_model.Maximize(x + 2 * y); + + // Set parameters, e.g. turn on logging. + SolveArguments args; + args.parameters.enable_output = true; + + // Solve and ensure an optimal solution was found with no errors. + const absl::StatusOr result = + Solve(lp_model, SolverType::kXpress, args); + CHECK_OK(result.status()); + CHECK_OK(result->termination.EnsureIsOptimal()); + + EXPECT_EQ(result->objective_value(), 2.5); + EXPECT_EQ(result->variable_values().at(x), 0.5); + EXPECT_EQ(result->variable_values().at(y), 1); +} + +SimpleLpTestParameters XpressDefaults() { + return SimpleLpTestParameters( + SolverType::kXpress, SolveParameters(), /*supports_duals=*/true, + /*supports_basis=*/true, + /*ensures_primal_ray=*/false, /*ensures_dual_ray=*/false, + /*disallows_infeasible_or_unbounded=*/true); +} + +INSTANTIATE_TEST_SUITE_P(XpressSolverLpTest, SimpleLpTest, + testing::Values(XpressDefaults())); + +int main(int argc, char** argv) { + testing::InitGoogleTest(&argc, argv); + RUN_ALL_TESTS(); +} \ No newline at end of file diff --git a/ortools/xpress/environment.h b/ortools/xpress/environment.h index 8bcb3a10767..ab61e07fcad 100644 --- a/ortools/xpress/environment.h +++ b/ortools/xpress/environment.h @@ -433,6 +433,26 @@ absl::Status LoadXpressDynamicLibrary(std::string& xpresspath); #define XPRS_MIP_UNBOUNDED 7 #define XPRS_OBJ_MINIMIZE 1 #define XPRS_OBJ_MAXIMIZE -1 +// *************************************************************************** +// * variable types * +// *************************************************************************** +#define XPRS_INTEGER 'I' +#define XPRS_CONTINUOUS 'C' +// *************************************************************************** +// * constraint types * +// *************************************************************************** +#define XPRS_LESS_EQUAL 'L' +#define XPRS_GREATER_EQUAL 'G' +#define XPRS_EQUAL 'E' +#define XPRS_RANGE 'R' +#define XPRS_NONBINDING 'N' +// *************************************************************************** +// * basis status * +// *************************************************************************** +#define XPRS_AT_LOWER 0 +#define XPRS_BASIC 1 +#define XPRS_AT_UPPER 2 +#define XPRS_FREE_SUPER 3 // Let's not reformat for rest of the file. // clang-format off From 073e8e3eaa2b85d01d6785b1267143b51aba444e Mon Sep 17 00:00:00 2001 From: Peter Mitri Date: Tue, 3 Dec 2024 16:19:22 +0100 Subject: [PATCH 02/14] changes after review --- ortools/math_opt/solvers/CMakeLists.txt | 6 +- ortools/math_opt/solvers/xpress/g_xpress.cc | 72 +++++++++---------- ortools/math_opt/solvers/xpress/g_xpress.h | 19 ++--- ortools/math_opt/solvers/xpress_solver.cc | 69 ++++++++++++------ ortools/math_opt/solvers/xpress_solver.h | 20 ++---- .../math_opt/solvers/xpress_solver_test.cc | 41 +---------- ortools/xpress/environment.h | 6 ++ 7 files changed, 107 insertions(+), 126 deletions(-) diff --git a/ortools/math_opt/solvers/CMakeLists.txt b/ortools/math_opt/solvers/CMakeLists.txt index a10c0ad3cc6..f6af550c765 100644 --- a/ortools/math_opt/solvers/CMakeLists.txt +++ b/ortools/math_opt/solvers/CMakeLists.txt @@ -40,6 +40,11 @@ if(NOT USE_SCIP) list(FILTER _SRCS EXCLUDE REGEX "/gscip_.*.h$") list(FILTER _SRCS EXCLUDE REGEX "/gscip_.*.cc$") endif() +if(NOT USE_XPRESS) + list(FILTER _SRCS EXCLUDE REGEX "/xpress/") + list(FILTER _SRCS EXCLUDE REGEX "/xpress_.*.h$") + list(FILTER _SRCS EXCLUDE REGEX "/xpress_.*.cc$") +endif() target_sources(${NAME} PRIVATE ${_SRCS}) set_target_properties(${NAME} PROPERTIES POSITION_INDEPENDENT_CODE ON) target_include_directories(${NAME} PUBLIC @@ -235,7 +240,6 @@ if(USE_HIGHS) endif() if(USE_XPRESS) - add_subdirectory(xpress) ortools_cxx_test( NAME math_opt_solvers_xpress_solver_test diff --git a/ortools/math_opt/solvers/xpress/g_xpress.cc b/ortools/math_opt/solvers/xpress/g_xpress.cc index ee6c8680585..efb060f4a2e 100644 --- a/ortools/math_opt/solvers/xpress/g_xpress.cc +++ b/ortools/math_opt/solvers/xpress/g_xpress.cc @@ -32,15 +32,6 @@ #include "ortools/base/status_macros.h" #include "ortools/xpress/environment.h" -// The argument to this macro is the invocation of a XPRS function that -// returns a status. If the function returns non-zero the macro aborts -// the program with an appropriate error message. -#define CHECK_STATUS(s) \ - do { \ - int const status_ = s; \ - CHECK_EQ(0, status_); \ - } while (0) - namespace operations_research::math_opt { constexpr int kXpressOk = 0; @@ -55,28 +46,31 @@ absl::Status Xpress::ToStatus(const int xprs_err, << "Xpress error code: " << xprs_err << ", message: " << errmsg; } -Xpress::Xpress(XPRSprob& model) - : xpress_model_(ABSL_DIE_IF_NULL(std::move(model))) {} +Xpress::Xpress(XPRSprob& model) : xpress_model_(ABSL_DIE_IF_NULL(model)) {} absl::StatusOr> Xpress::New( const std::string& model_name) { bool correctlyLoaded = initXpressEnv(); CHECK(correctlyLoaded); XPRSprob* model; - CHECK_STATUS(XPRScreateprob(model)); + CHECK_EQ(kXpressOk, XPRScreateprob(model)); DCHECK(model != nullptr); // should not be NULL if status=0 - CHECK_STATUS(XPRSaddcbmessage(*model, optimizermsg, NULL, 0)); + CHECK_EQ(kXpressOk, XPRSaddcbmessage(*model, printXpressMessage, NULL, 0)); return absl::WrapUnique(new Xpress(*model)); } -void XPRS_CC optimizermsg(XPRSprob prob, void* data, const char* sMsg, int nLen, - int nMsgLvl) { - printf("%*s\n", nLen, sMsg); +void XPRS_CC Xpress::printXpressMessage(XPRSprob prob, void* data, + const char* sMsg, int nLen, + int nMsgLvl) { + if (sMsg) { + const std::string strMsg = sMsg; + std::cout << strMsg << std::endl; + } } Xpress::~Xpress() { - CHECK_STATUS(XPRSdestroyprob(xpress_model_)); - CHECK_STATUS(XPRSfree()); + CHECK_EQ(kXpressOk, XPRSdestroyprob(xpress_model_)); + CHECK_EQ(kXpressOk, XPRSfree()); } absl::Status Xpress::AddVars(const absl::Span obj, @@ -105,9 +99,9 @@ absl::Status Xpress::AddVars(const absl::Span vbegin, if (!vbegin.empty()) { CHECK_EQ(vbegin.size(), num_vars); } - - CHECK_STATUS(XPRSaddcols(xpress_model_, num_vars, 0, c_obj, nullptr, nullptr, - nullptr, lb.data(), ub.data())); + // TODO: look into int64 support for number of vars (use XPRSaddcols64) + CHECK_EQ(kXpressOk, XPRSaddcols(xpress_model_, num_vars, 0, c_obj, nullptr, + nullptr, nullptr, lb.data(), ub.data())); return absl::OkStatus(); } @@ -122,27 +116,28 @@ absl::Status Xpress::AddConstrs(const absl::Span sense, absl::Status Xpress::SetObjective(bool maximize, double offset, const absl::Span colind, const absl::Span values) { - int n_cols = colind.size(); + const int n_cols = static_cast(colind.size()); auto c_colind = const_cast(colind.data()); auto c_values = const_cast(values.data()); - CHECK_STATUS(XPRSchgobj(xpress_model_, n_cols, c_colind, c_values)); + CHECK_EQ(kXpressOk, XPRSchgobj(xpress_model_, n_cols, c_colind, c_values)); static int indexes[1] = {-1}; double xprs_values[1] = {-offset}; - CHECK_STATUS(XPRSchgobj(xpress_model_, 1, indexes, xprs_values)); - CHECK_STATUS(XPRSchgobjsense( - xpress_model_, maximize ? XPRS_OBJ_MAXIMIZE : XPRS_OBJ_MINIMIZE)); + CHECK_EQ(kXpressOk, XPRSchgobj(xpress_model_, 1, indexes, xprs_values)); + CHECK_EQ(kXpressOk, + XPRSchgobjsense(xpress_model_, + maximize ? XPRS_OBJ_MAXIMIZE : XPRS_OBJ_MINIMIZE)); return absl::OkStatus(); } absl::Status Xpress::ChgCoeffs(absl::Span rowind, absl::Span colind, absl::Span values) { - int n_coefs = rowind.size(); + const int n_coefs = static_cast(rowind.size()); auto c_rowind = const_cast(rowind.data()); auto c_colind = const_cast(colind.data()); auto c_values = const_cast(values.data()); - CHECK_STATUS( - XPRSchgmcoef(xpress_model_, n_coefs, c_rowind, c_colind, c_values)); + CHECK_EQ(kXpressOk, + XPRSchgmcoef(xpress_model_, n_coefs, c_rowind, c_colind, c_values)); return absl::OkStatus(); } @@ -192,12 +187,11 @@ absl::StatusOr Xpress::GetDoubleAttr(int attribute) const { absl::StatusOr> Xpress::GetPrimalValues() const { int nCols = GetNumberOfColumns(); XPRSgetintattrib(xpress_model_, XPRS_COLS, &nCols); - double values[nCols]; - RETURN_IF_ERROR( - ToStatus(XPRSgetlpsol(xpress_model_, values, nullptr, nullptr, nullptr))) + std::vector values(nCols); + RETURN_IF_ERROR(ToStatus( + XPRSgetlpsol(xpress_model_, values.data(), nullptr, nullptr, nullptr))) << "Error getting Xpress LP solution"; - std::vector result(values, values + nCols); - return std::move(result); + return values; } int Xpress::GetNumberOfRows() const { @@ -215,21 +209,23 @@ int Xpress::GetNumberOfColumns() const { std::vector Xpress::GetConstraintDuals() const { int nRows = GetNumberOfRows(); double values[nRows]; - CHECK_STATUS(XPRSgetlpsol(xpress_model_, nullptr, nullptr, values, nullptr)); + CHECK_EQ(kXpressOk, + XPRSgetlpsol(xpress_model_, nullptr, nullptr, values, nullptr)); return {values, values + nRows}; } std::vector Xpress::GetReducedCostValues() const { int nCols = GetNumberOfColumns(); double values[nCols]; - CHECK_STATUS(XPRSgetlpsol(xpress_model_, nullptr, nullptr, nullptr, values)); + CHECK_EQ(kXpressOk, + XPRSgetlpsol(xpress_model_, nullptr, nullptr, nullptr, values)); return {values, values + nCols}; } std::vector Xpress::GetVariableBasis() const { int const nCols = GetNumberOfColumns(); std::vector basis(nCols); - CHECK_STATUS(XPRSgetbasis(xpress_model_, basis.data(), 0)); + CHECK_EQ(kXpressOk, XPRSgetbasis(xpress_model_, basis.data(), 0)); return basis; } -} // namespace operations_research::math_opt \ No newline at end of file +} // namespace operations_research::math_opt diff --git a/ortools/math_opt/solvers/xpress/g_xpress.h b/ortools/math_opt/solvers/xpress/g_xpress.h index 9d5fee2909b..0e136b3d978 100644 --- a/ortools/math_opt/solvers/xpress/g_xpress.h +++ b/ortools/math_opt/solvers/xpress/g_xpress.h @@ -42,10 +42,6 @@ class Xpress { public: Xpress() = delete; - absl::Status ToStatus( - int xprs_err, - absl::StatusCode code = absl::StatusCode::kInvalidArgument) const; - // Creates a new Xpress static absl::StatusOr> New( const std::string& model_name); @@ -87,16 +83,21 @@ class Xpress { std::vector GetReducedCostValues() const; std::vector GetVariableBasis() const; + static void XPRS_CC printXpressMessage(XPRSprob prob, void* data, + const char* sMsg, int nLen, + int nMsgLvl); + + int GetNumberOfRows() const; + int GetNumberOfColumns() const; + private: XPRSprob xpress_model_; explicit Xpress(XPRSprob& model); - int GetNumberOfRows() const; - int GetNumberOfColumns() const; - - static void XPRS_CC optimizermsg(XPRSprob prob, void* data, const char* sMsg, - int nLen, int nMsgLvl); + absl::Status ToStatus( + int xprs_err, + absl::StatusCode code = absl::StatusCode::kInvalidArgument) const; }; } // namespace operations_research::math_opt diff --git a/ortools/math_opt/solvers/xpress_solver.cc b/ortools/math_opt/solvers/xpress_solver.cc index 71e92d72d93..2cd717ea668 100644 --- a/ortools/math_opt/solvers/xpress_solver.cc +++ b/ortools/math_opt/solvers/xpress_solver.cc @@ -11,7 +11,7 @@ // See the License for the specific language governing permissions and // limitations under the License. -#include "xpress_solver.h" +#include "ortools/math_opt/solvers/xpress_solver.h" #include "ortools/base/map_util.h" #include "ortools/base/protoutil.h" @@ -66,9 +66,10 @@ absl::Status XpressSolver::AddNewVariables( const VariablesProto& new_variables) { const int num_new_variables = new_variables.lower_bounds().size(); std::vector variable_type(num_new_variables); + int n_variables = xpress_->GetNumberOfColumns(); for (int j = 0; j < num_new_variables; ++j) { - const VariableId id = new_variables.ids(j); - InsertOrDie(&variables_map_, id, j + num_xpress_variables_); + const VarId id = new_variables.ids(j); + InsertOrDie(&variables_map_, id, j + n_variables); variable_type[j] = new_variables.integers(j) ? XPRS_INTEGER : XPRS_CONTINUOUS; if (new_variables.integers(j)) { @@ -79,7 +80,6 @@ absl::Status XpressSolver::AddNewVariables( RETURN_IF_ERROR(xpress_->AddVars({}, new_variables.lower_bounds(), new_variables.upper_bounds(), variable_type)); - num_xpress_variables_ += num_new_variables; // Not adding names for performance (have to call XPRSaddnames) // TODO : keep names in a cache and add them when needed? @@ -109,6 +109,7 @@ absl::Status XpressSolver::AddNewLinearConstraints( constraint_rhs.reserve(num_new_constraints); constraint_sense.reserve(num_new_constraints); new_slacks.reserve(num_new_constraints); + int n_constraints = xpress_->GetNumberOfRows(); for (int i = 0; i < num_new_constraints; ++i) { const int64_t id = constraints.ids(i); LinearConstraintData& constraint_data = @@ -117,7 +118,7 @@ absl::Status XpressSolver::AddNewLinearConstraints( const double ub = constraints.upper_bounds(i); constraint_data.lower_bound = lb; constraint_data.upper_bound = ub; - constraint_data.constraint_index = i + num_xpress_lin_cons_; + constraint_data.constraint_index = i + n_constraints; char sense = XPRS_EQUAL; double rhs = 0.0; const bool lb_is_xprs_neg_inf = lb <= XPRS_MINUSINFINITY; @@ -134,7 +135,7 @@ absl::Status XpressSolver::AddNewLinearConstraints( } else { // Note that constraints where the lower bound and the upper bound are // -+infinity translate into a range constraint with an unbounded slack. - constraint_data.slack_index = new_slacks.size() + num_xpress_variables_; + constraint_data.slack_index = new_slacks.size() + n_constraints; new_slacks.push_back(&constraint_data); } constraint_rhs.emplace_back(rhs); @@ -142,7 +143,6 @@ absl::Status XpressSolver::AddNewLinearConstraints( } // Add all constraints in one call. RETURN_IF_ERROR(xpress_->AddConstrs(constraint_sense, constraint_rhs)); - num_xpress_lin_cons_ += num_new_constraints; // Add slacks for true ranged constraints (if needed) if (!new_slacks.empty()) { // TODO @@ -153,6 +153,8 @@ absl::Status XpressSolver::AddNewLinearConstraints( } absl::Status XpressSolver::AddSingleObjective(const ObjectiveProto& objective) { + // TODO: reactivate the following code after figuring out why the condition is + // true for LP tests /*if (objective.has_quadratic_coefficients()) { return absl::UnimplementedError( "Quadratic objectives are not yet implemented in XPRESS solver " @@ -176,12 +178,14 @@ absl::Status XpressSolver::AddSingleObjective(const ObjectiveProto& objective) { absl::Status XpressSolver::ChangeCoefficients( const SparseDoubleMatrixProto& matrix) { const int num_coefficients = matrix.row_ids().size(); - std::vector row_index(num_coefficients); - std::vector col_index(num_coefficients); + std::vector row_index; + row_index.reserve(num_coefficients); + std::vector col_index; + col_index.reserve(num_coefficients); for (int k = 0; k < num_coefficients; ++k) { - row_index[k] = - linear_constraints_map_.at(matrix.row_ids(k)).constraint_index; - col_index[k] = variables_map_.at(matrix.column_ids(k)); + row_index.push_back( + linear_constraints_map_.at(matrix.row_ids(k)).constraint_index); + col_index.push_back(variables_map_.at(matrix.column_ids(k))); } return xpress_->ChgCoeffs(row_index, col_index, matrix.coefficients()); } @@ -250,20 +254,15 @@ absl::StatusOr XpressSolver::ExtractSolveResultProto( ASSIGN_OR_RETURN(*result.mutable_termination(), ConvertTerminationReason(solution_claims, best_primal_bound, best_dual_bound)); - return std::move(result); + return result; } absl::StatusOr XpressSolver::GetBestPrimalBound() const { - ASSIGN_OR_RETURN( - const double best_primal_bound, - xpress_->GetDoubleAttr(is_mip_ ? XPRS_MIPOBJVAL : XPRS_LPOBJVAL)); - return best_primal_bound; + return xpress_->GetDoubleAttr(is_mip_ ? XPRS_MIPOBJVAL : XPRS_LPOBJVAL); } absl::StatusOr XpressSolver::GetBestDualBound() const { - ASSIGN_OR_RETURN(const double best_dual_bound, - xpress_->GetDoubleAttr(XPRS_BESTBOUND)); - return best_dual_bound; + return xpress_->GetDoubleAttr(XPRS_BESTBOUND); } absl::StatusOr XpressSolver::GetSolutions( @@ -385,7 +384,7 @@ XpressSolver::GetConvexDualSolutionIfAvailable( "unavailable or infinite, and no dual feasible solution is returned"); } return SolutionAndClaim{ - .solution = std::move(dual_solution), + .solution = dual_solution, .feasible_solution_exists = dual_feasible_solution_exists}; } @@ -432,7 +431,7 @@ absl::StatusOr> XpressSolver::GetBasisIfAvailable() { } else if (xpress_status_ == XPRS_LP_UNBOUNDED) { basis.set_basic_dual_feasibility(SOLUTION_STATUS_INFEASIBLE); } - return std::move(basis); + return basis; } absl::StatusOr XpressSolver::GetSolveStats( @@ -465,6 +464,10 @@ absl::StatusOr XpressSolver::ConvertTerminationReason( // TODO : improve this if (!is_mip_) { switch (xpress_status_) { + case XPRS_LP_UNSTARTED: + return TerminateForReason( + is_maximize_, TERMINATION_REASON_OTHER_ERROR, + "Problem solve has not started (XPRS_LP_UNSTARTED)"); case XPRS_LP_OPTIMAL: return OptimalTerminationProto(best_primal_bound, best_dual_bound); case XPRS_LP_INFEAS: @@ -472,6 +475,13 @@ absl::StatusOr XpressSolver::ConvertTerminationReason( is_maximize_, solution_claims.dual_feasible_solution_exists ? FEASIBILITY_STATUS_FEASIBLE : FEASIBILITY_STATUS_UNDETERMINED); + case XPRS_LP_CUTOFF: + return CutoffTerminationProto( + is_maximize_, "Objective worse than cutoff (XPRS_LP_CUTOFF)"); + case XPRS_LP_UNFINISHED: + return LimitTerminationProto( + is_maximize_, LIMIT_UNSPECIFIED, best_primal_bound, best_dual_bound, + "Solve did not finish (XPRS_LP_UNFINISHED)"); case XPRS_LP_UNBOUNDED: if (solution_claims.primal_feasible_solution_exists) { return UnboundedTerminationProto(is_maximize_); @@ -479,8 +489,21 @@ absl::StatusOr XpressSolver::ConvertTerminationReason( return InfeasibleOrUnboundedTerminationProto( is_maximize_, FEASIBILITY_STATUS_INFEASIBLE, "Xpress status XPRS_LP_UNBOUNDED"); + case XPRS_LP_CUTOFF_IN_DUAL: + return CutoffTerminationProto( + is_maximize_, "Cutoff in dual (XPRS_LP_CUTOFF_IN_DUAL)"); + case XPRS_LP_UNSOLVED: + return TerminateForReason( + is_maximize_, TERMINATION_REASON_NUMERICAL_ERROR, + "Problem could not be solved due to numerical issues " + "(XPRS_LP_UNSOLVED)"); + case XPRS_LP_NONCONVEX: + return TerminateForReason(is_maximize_, TERMINATION_REASON_OTHER_ERROR, + "Problem contains quadratic data, which is " + "not convex (XPRS_LP_NONCONVEX)"); default: - return TerminateForReason(is_maximize_, TERMINATION_REASON_IMPRECISE); + return absl::InternalError(absl::StrCat( + "Missing Xpress LP status code case: ", xpress_status_)); } } else { return absl::UnimplementedError("XpressSolver does not handle MIPs yet"); diff --git a/ortools/math_opt/solvers/xpress_solver.h b/ortools/math_opt/solvers/xpress_solver.h index 4b250c09f51..5474d829da5 100644 --- a/ortools/math_opt/solvers/xpress_solver.h +++ b/ortools/math_opt/solvers/xpress_solver.h @@ -71,7 +71,7 @@ class XpressSolver : public SolverInterface { explicit XpressSolver(std::unique_ptr g_xpress); // For easing reading the code, we declare these types: - using VariableId = int64_t; + using VarId = int64_t; using AuxiliaryObjectiveId = int64_t; using LinearConstraintId = int64_t; using QuadraticConstraintId = int64_t; @@ -111,11 +111,12 @@ class XpressSolver : public SolverInterface { }; struct SolutionClaims { - bool primal_feasible_solution_exists; - bool dual_feasible_solution_exists; + bool primal_feasible_solution_exists = false; + bool dual_feasible_solution_exists = false; }; struct SolutionsAndClaims { + // TODO: simplify this structure or change it std::vector solutions; SolutionClaims solution_claims; }; @@ -181,7 +182,7 @@ class XpressSolver : public SolverInterface { // Internal correspondence from variable proto IDs to Xpress-numbered // variables. - gtl::linked_hash_map variables_map_; + gtl::linked_hash_map variables_map_; // Internal correspondence from linear constraint proto IDs to // Xpress-numbered linear constraint and extra information. gtl::linked_hash_map @@ -198,17 +199,6 @@ class XpressSolver : public SolverInterface { // quantities are updated immediately after adding or removing to the model, // so it is correct even if XPRESS C API has not yet been called. - // Number of Xpress variables. - int num_xpress_variables_ = 0; - // Number of Xpress linear constraints. - int num_xpress_lin_cons_ = 0; - // Number of Xpress quadratic constraints. - int num_xpress_quad_cons_ = 0; - // Number of Xpress SOS constraints. - int num_xpress_sos_cons_ = 0; - // Number of Xpress general constraints. - int num_xpress_gen_cons_ = 0; - bool is_mip_ = false; bool is_maximize_ = false; int xpress_status_ = 0; diff --git a/ortools/math_opt/solvers/xpress_solver_test.cc b/ortools/math_opt/solvers/xpress_solver_test.cc index 87eadab6ab5..107bd9e940a 100644 --- a/ortools/math_opt/solvers/xpress_solver_test.cc +++ b/ortools/math_opt/solvers/xpress_solver_test.cc @@ -29,45 +29,6 @@ using namespace operations_research::math_opt; -using ::testing::AnyOf; -using ::testing::status::IsOkAndHolds; -constexpr double kTolerance = 1.0e-5; -constexpr double kInf = std::numeric_limits::infinity(); - -struct SolvedModel { - std::unique_ptr model; - SolveResult expected_result; -}; - -// max(x + 2 * y) -// s.t.: -1 <= x <= 1.5 -// 0 <= y <= 1 -// x + y <= 1.5 -// optimal solution is: (x, y) = (0.5, 1) -// obj_value = 2.5 -TEST(XpressSolverTest, SimpleLpTwoVar) { - // Build the model. - Model lp_model("getting_started_lp"); - const Variable x = lp_model.AddContinuousVariable(-1.0, 1.5, "x"); - const Variable y = lp_model.AddContinuousVariable(0.0, 1.0, "y"); - lp_model.AddLinearConstraint(x + y <= 1.5, "c"); - lp_model.Maximize(x + 2 * y); - - // Set parameters, e.g. turn on logging. - SolveArguments args; - args.parameters.enable_output = true; - - // Solve and ensure an optimal solution was found with no errors. - const absl::StatusOr result = - Solve(lp_model, SolverType::kXpress, args); - CHECK_OK(result.status()); - CHECK_OK(result->termination.EnsureIsOptimal()); - - EXPECT_EQ(result->objective_value(), 2.5); - EXPECT_EQ(result->variable_values().at(x), 0.5); - EXPECT_EQ(result->variable_values().at(y), 1); -} - SimpleLpTestParameters XpressDefaults() { return SimpleLpTestParameters( SolverType::kXpress, SolveParameters(), /*supports_duals=*/true, @@ -82,4 +43,4 @@ INSTANTIATE_TEST_SUITE_P(XpressSolverLpTest, SimpleLpTest, int main(int argc, char** argv) { testing::InitGoogleTest(&argc, argv); RUN_ALL_TESTS(); -} \ No newline at end of file +} diff --git a/ortools/xpress/environment.h b/ortools/xpress/environment.h index ab61e07fcad..ae6e34a7c1a 100644 --- a/ortools/xpress/environment.h +++ b/ortools/xpress/environment.h @@ -424,9 +424,15 @@ absl::Status LoadXpressDynamicLibrary(std::string& xpresspath); #define XPRS_MIPSTATUS 1011 #define XPRS_NODES 1013 #define XPRS_COLS 1018 +#define XPRS_LP_UNSTARTED 0 #define XPRS_LP_OPTIMAL 1 #define XPRS_LP_INFEAS 2 +#define XPRS_LP_CUTOFF 3 +#define XPRS_LP_UNFINISHED 4 #define XPRS_LP_UNBOUNDED 5 +#define XPRS_LP_CUTOFF_IN_DUAL 6 +#define XPRS_LP_UNSOLVED 7 +#define XPRS_LP_NONCONVEX 8 #define XPRS_MIP_SOLUTION 4 #define XPRS_MIP_INFEAS 5 #define XPRS_MIP_OPTIMAL 6 From 6a6d3c094b835fdebea8866d9f41c17563d83391 Mon Sep 17 00:00:00 2001 From: Peter Mitri Date: Tue, 3 Dec 2024 16:21:26 +0100 Subject: [PATCH 03/14] simplification --- ortools/math_opt/solvers/xpress/g_xpress.cc | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/ortools/math_opt/solvers/xpress/g_xpress.cc b/ortools/math_opt/solvers/xpress/g_xpress.cc index efb060f4a2e..dfcafa9ec65 100644 --- a/ortools/math_opt/solvers/xpress/g_xpress.cc +++ b/ortools/math_opt/solvers/xpress/g_xpress.cc @@ -63,8 +63,7 @@ void XPRS_CC Xpress::printXpressMessage(XPRSprob prob, void* data, const char* sMsg, int nLen, int nMsgLvl) { if (sMsg) { - const std::string strMsg = sMsg; - std::cout << strMsg << std::endl; + std::cout << sMsg << std::endl; } } From 95357b33eac53f0efe4a60eb4e78e731dfd932f6 Mon Sep 17 00:00:00 2001 From: Peter Mitri Date: Tue, 3 Dec 2024 17:21:46 +0100 Subject: [PATCH 04/14] add support for ranged constraints --- ortools/math_opt/solvers/xpress/g_xpress.cc | 5 +-- ortools/math_opt/solvers/xpress/g_xpress.h | 4 ++- ortools/math_opt/solvers/xpress_solver.cc | 37 ++++++++------------- 3 files changed, 19 insertions(+), 27 deletions(-) diff --git a/ortools/math_opt/solvers/xpress/g_xpress.cc b/ortools/math_opt/solvers/xpress/g_xpress.cc index dfcafa9ec65..4c00b3de073 100644 --- a/ortools/math_opt/solvers/xpress/g_xpress.cc +++ b/ortools/math_opt/solvers/xpress/g_xpress.cc @@ -105,11 +105,12 @@ absl::Status Xpress::AddVars(const absl::Span vbegin, } absl::Status Xpress::AddConstrs(const absl::Span sense, - const absl::Span rhs) { + const absl::Span rhs, + const absl::Span rng) { const int num_cons = static_cast(sense.size()); CHECK_EQ(rhs.size(), num_cons); return ToStatus(XPRSaddrows(xpress_model_, num_cons, 0, sense.data(), - rhs.data(), NULL, NULL, NULL, NULL)); + rhs.data(), rng.data(), NULL, NULL, NULL)); } absl::Status Xpress::SetObjective(bool maximize, double offset, diff --git a/ortools/math_opt/solvers/xpress/g_xpress.h b/ortools/math_opt/solvers/xpress/g_xpress.h index 0e136b3d978..f5adbcd676b 100644 --- a/ortools/math_opt/solvers/xpress/g_xpress.h +++ b/ortools/math_opt/solvers/xpress/g_xpress.h @@ -64,7 +64,9 @@ class Xpress { absl::Span vtype); absl::Status AddConstrs(absl::Span sense, - absl::Span rhs); + absl::Span rhs, + absl::Span rng); + absl::Status SetObjective(bool maximize, double offset, absl::Span colind, absl::Span values); diff --git a/ortools/math_opt/solvers/xpress_solver.cc b/ortools/math_opt/solvers/xpress_solver.cc index 2cd717ea668..f70ef891758 100644 --- a/ortools/math_opt/solvers/xpress_solver.cc +++ b/ortools/math_opt/solvers/xpress_solver.cc @@ -94,17 +94,9 @@ absl::Status XpressSolver::AddNewLinearConstraints( const LinearConstraintsProto& constraints) { const int num_new_constraints = constraints.lower_bounds().size(); - // Constraints are translated into: - // 1. ax <= upper_bound (if lower bound <= -INFINITY, and upper_bound - // is finite and less than INFINITY) - // 2. ax >= lower_bound (if upper bound >= INFINITY, and lower_bound is - // finite and greater than -INFINITY) - // 3. ax == xxxxx_bound (if both bounds are finite, equal, and their - // absolute values less than INFINITY) - // 4. ax - slack = 0.0 (otherwise, - // slack bounds == [lower_bound, upper_bound]) - std::vector constraint_rhs; std::vector constraint_sense; + std::vector constraint_rhs; + std::vector constraint_rng; std::vector new_slacks; constraint_rhs.reserve(num_new_constraints); constraint_sense.reserve(num_new_constraints); @@ -113,7 +105,7 @@ absl::Status XpressSolver::AddNewLinearConstraints( for (int i = 0; i < num_new_constraints; ++i) { const int64_t id = constraints.ids(i); LinearConstraintData& constraint_data = - gtl::InsertKeyOrDie(&linear_constraints_map_, id); + InsertKeyOrDie(&linear_constraints_map_, id); const double lb = constraints.lower_bounds(i); const double ub = constraints.upper_bounds(i); constraint_data.lower_bound = lb; @@ -121,6 +113,7 @@ absl::Status XpressSolver::AddNewLinearConstraints( constraint_data.constraint_index = i + n_constraints; char sense = XPRS_EQUAL; double rhs = 0.0; + double rng = 0.0; const bool lb_is_xprs_neg_inf = lb <= XPRS_MINUSINFINITY; const bool ub_is_xprs_pos_inf = ub >= XPRS_PLUSINFINITY; if (lb_is_xprs_neg_inf && !ub_is_xprs_pos_inf) { @@ -133,23 +126,19 @@ absl::Status XpressSolver::AddNewLinearConstraints( sense = XPRS_EQUAL; rhs = lb; } else { - // Note that constraints where the lower bound and the upper bound are - // -+infinity translate into a range constraint with an unbounded slack. - constraint_data.slack_index = new_slacks.size() + n_constraints; - new_slacks.push_back(&constraint_data); + if (ub < lb) { + return absl::InvalidArgumentError("Lower bound > Upper bound"); + } + sense = XPRS_RANGE; + rhs = ub; + rng = ub - lb; } - constraint_rhs.emplace_back(rhs); constraint_sense.emplace_back(sense); + constraint_rhs.emplace_back(rhs); + constraint_rng.emplace_back(rng); } // Add all constraints in one call. - RETURN_IF_ERROR(xpress_->AddConstrs(constraint_sense, constraint_rhs)); - // Add slacks for true ranged constraints (if needed) - if (!new_slacks.empty()) { - // TODO - return absl::UnimplementedError( - "XpressSolver does not support slack variables yet"); - } - return absl::OkStatus(); + return xpress_->AddConstrs(constraint_sense, constraint_rhs, constraint_rng); } absl::Status XpressSolver::AddSingleObjective(const ObjectiveProto& objective) { From 19805c0bc5c3905538d447141bf839a53cd5018a Mon Sep 17 00:00:00 2001 From: Peter Mitri Date: Tue, 3 Dec 2024 18:08:37 +0100 Subject: [PATCH 05/14] fix some tests --- ortools/math_opt/solvers/xpress_solver.cc | 18 +++++++++++------- 1 file changed, 11 insertions(+), 7 deletions(-) diff --git a/ortools/math_opt/solvers/xpress_solver.cc b/ortools/math_opt/solvers/xpress_solver.cc index f70ef891758..78ec2ee6374 100644 --- a/ortools/math_opt/solvers/xpress_solver.cc +++ b/ortools/math_opt/solvers/xpress_solver.cc @@ -251,7 +251,8 @@ absl::StatusOr XpressSolver::GetBestPrimalBound() const { } absl::StatusOr XpressSolver::GetBestDualBound() const { - return xpress_->GetDoubleAttr(XPRS_BESTBOUND); + // TODO: setting LP primal value as best dual bound. Can this be improved? + return xpress_->GetDoubleAttr(XPRS_LPOBJVAL); } absl::StatusOr XpressSolver::GetSolutions( @@ -302,6 +303,7 @@ bool XpressSolver::isFeasible() const { } SolutionStatusProto XpressSolver::getLpSolutionStatus() const { + // TODO : put all statuses here switch (xpress_status_) { case XPRS_LP_OPTIMAL: return SOLUTION_STATUS_FEASIBLE; @@ -318,18 +320,18 @@ absl::StatusOr> XpressSolver::GetConvexPrimalSolutionIfAvailable( const ModelSolveParametersProto& model_parameters) const { PrimalSolutionProto primal_solution; + primal_solution.set_feasibility_status(getLpSolutionStatus()); if (isFeasible()) { ASSIGN_OR_RETURN(const double sol_val, xpress_->GetDoubleAttr(XPRS_LPOBJVAL)); primal_solution.set_objective_value(sol_val); + XpressVectorToSparseDoubleVector(xpress_->GetPrimalValues().value(), + variables_map_, + *primal_solution.mutable_variable_values(), + model_parameters.variable_values_filter()); } else { // TODO } - primal_solution.set_feasibility_status(getLpSolutionStatus()); - XpressVectorToSparseDoubleVector(xpress_->GetPrimalValues().value(), - variables_map_, - *primal_solution.mutable_variable_values(), - model_parameters.variable_values_filter()); const bool primal_feasible_solution_exists = (primal_solution.feasibility_status() == SOLUTION_STATUS_FEASIBLE); return SolutionAndClaim{ @@ -364,7 +366,9 @@ XpressSolver::GetConvexDualSolutionIfAvailable( dual_solution.set_feasibility_status(getLpSolutionStatus()); bool dual_feasible_solution_exists = (dual_solution.feasibility_status() == SOLUTION_STATUS_FEASIBLE); - ASSIGN_OR_RETURN(const double best_dual_bound, GetBestDualBound()); + ASSIGN_OR_RETURN(const double best_dual_bound, + xpress_->GetDoubleAttr(XPRS_LPOBJVAL)); + // const double best_dual_bound = is_maximize_ ? kMinusInf : kPlusInf; if (dual_feasible_solution_exists || std::isfinite(best_dual_bound)) { dual_feasible_solution_exists = true; } else if (xpress_status_ == XPRS_LP_OPTIMAL) { From d15a1b6509ea9d61a463ee34fb3171495d0426a3 Mon Sep 17 00:00:00 2001 From: Peter Mitri Date: Wed, 4 Dec 2024 11:14:17 +0100 Subject: [PATCH 06/14] Fix basis handling, and infeasible & unbounded LPs handling. Refactor some code. --- ortools/math_opt/solvers/xpress/g_xpress.cc | 61 ++++++++------- ortools/math_opt/solvers/xpress/g_xpress.h | 6 +- ortools/math_opt/solvers/xpress_solver.cc | 76 ++++++++++--------- .../math_opt/solvers/xpress_solver_test.cc | 36 ++++++++- 4 files changed, 114 insertions(+), 65 deletions(-) diff --git a/ortools/math_opt/solvers/xpress/g_xpress.cc b/ortools/math_opt/solvers/xpress/g_xpress.cc index 4c00b3de073..838a69eb878 100644 --- a/ortools/math_opt/solvers/xpress/g_xpress.cc +++ b/ortools/math_opt/solvers/xpress/g_xpress.cc @@ -32,6 +32,8 @@ #include "ortools/base/status_macros.h" #include "ortools/xpress/environment.h" +// TODO : try to remove all CHECK_EQ and replace them with absl Status + namespace operations_research::math_opt { constexpr int kXpressOk = 0; @@ -99,9 +101,8 @@ absl::Status Xpress::AddVars(const absl::Span vbegin, CHECK_EQ(vbegin.size(), num_vars); } // TODO: look into int64 support for number of vars (use XPRSaddcols64) - CHECK_EQ(kXpressOk, XPRSaddcols(xpress_model_, num_vars, 0, c_obj, nullptr, - nullptr, nullptr, lb.data(), ub.data())); - return absl::OkStatus(); + return ToStatus(XPRSaddcols(xpress_model_, num_vars, 0, c_obj, nullptr, + nullptr, nullptr, lb.data(), ub.data())); } absl::Status Xpress::AddConstrs(const absl::Span sense, @@ -116,17 +117,19 @@ absl::Status Xpress::AddConstrs(const absl::Span sense, absl::Status Xpress::SetObjective(bool maximize, double offset, const absl::Span colind, const absl::Span values) { + RETURN_IF_ERROR(ToStatus(XPRSchgobjsense( + xpress_model_, maximize ? XPRS_OBJ_MAXIMIZE : XPRS_OBJ_MINIMIZE))) + << "Failed to change objective sense in XPRESS"; + + static int indexes[1] = {-1}; + double xprs_values[1] = {-offset}; + RETURN_IF_ERROR(ToStatus(XPRSchgobj(xpress_model_, 1, indexes, xprs_values))) + << "Failed to set objective offset in XPRESS"; + const int n_cols = static_cast(colind.size()); auto c_colind = const_cast(colind.data()); auto c_values = const_cast(values.data()); - CHECK_EQ(kXpressOk, XPRSchgobj(xpress_model_, n_cols, c_colind, c_values)); - static int indexes[1] = {-1}; - double xprs_values[1] = {-offset}; - CHECK_EQ(kXpressOk, XPRSchgobj(xpress_model_, 1, indexes, xprs_values)); - CHECK_EQ(kXpressOk, - XPRSchgobjsense(xpress_model_, - maximize ? XPRS_OBJ_MAXIMIZE : XPRS_OBJ_MINIMIZE)); - return absl::OkStatus(); + return ToStatus(XPRSchgobj(xpress_model_, n_cols, c_colind, c_values)); } absl::Status Xpress::ChgCoeffs(absl::Span rowind, @@ -136,9 +139,8 @@ absl::Status Xpress::ChgCoeffs(absl::Span rowind, auto c_rowind = const_cast(rowind.data()); auto c_colind = const_cast(colind.data()); auto c_values = const_cast(values.data()); - CHECK_EQ(kXpressOk, - XPRSchgmcoef(xpress_model_, n_coefs, c_rowind, c_colind, c_values)); - return absl::OkStatus(); + return ToStatus( + XPRSchgmcoef(xpress_model_, n_coefs, c_rowind, c_colind, c_values)); } absl::StatusOr Xpress::LpOptimizeAndGetStatus() { @@ -206,26 +208,31 @@ int Xpress::GetNumberOfColumns() const { return n; } -std::vector Xpress::GetConstraintDuals() const { +absl::StatusOr> Xpress::GetConstraintDuals() const { int nRows = GetNumberOfRows(); double values[nRows]; - CHECK_EQ(kXpressOk, - XPRSgetlpsol(xpress_model_, nullptr, nullptr, values, nullptr)); - return {values, values + nRows}; + RETURN_IF_ERROR( + ToStatus(XPRSgetlpsol(xpress_model_, nullptr, nullptr, values, nullptr))) + << "Failed to retrieve LP solution from XPRESS"; + std::vector result(values, values + nRows); + return result; } -std::vector Xpress::GetReducedCostValues() const { +absl::StatusOr> Xpress::GetReducedCostValues() const { int nCols = GetNumberOfColumns(); double values[nCols]; - CHECK_EQ(kXpressOk, - XPRSgetlpsol(xpress_model_, nullptr, nullptr, nullptr, values)); - return {values, values + nCols}; + RETURN_IF_ERROR( + ToStatus(XPRSgetlpsol(xpress_model_, nullptr, nullptr, nullptr, values))) + << "Failed to retrieve LP solution from XPRESS"; + std::vector result(values, values + nCols); + return result; } -std::vector Xpress::GetVariableBasis() const { - int const nCols = GetNumberOfColumns(); - std::vector basis(nCols); - CHECK_EQ(kXpressOk, XPRSgetbasis(xpress_model_, basis.data(), 0)); - return basis; +absl::Status Xpress::GetBasis(std::vector& rowBasis, + std::vector& colBasis) const { + rowBasis.resize(GetNumberOfRows()); + colBasis.resize(GetNumberOfColumns()); + return ToStatus( + XPRSgetbasis(xpress_model_, rowBasis.data(), colBasis.data())); } } // namespace operations_research::math_opt diff --git a/ortools/math_opt/solvers/xpress/g_xpress.h b/ortools/math_opt/solvers/xpress/g_xpress.h index f5adbcd676b..05419536256 100644 --- a/ortools/math_opt/solvers/xpress/g_xpress.h +++ b/ortools/math_opt/solvers/xpress/g_xpress.h @@ -81,9 +81,9 @@ class Xpress { void Terminate(); absl::StatusOr> GetPrimalValues() const; - std::vector GetConstraintDuals() const; - std::vector GetReducedCostValues() const; - std::vector GetVariableBasis() const; + absl::StatusOr> GetConstraintDuals() const; + absl::StatusOr> GetReducedCostValues() const; + absl::Status GetBasis(std::vector& rowBasis, std::vector& colBasis) const; static void XPRS_CC printXpressMessage(XPRSprob prob, void* data, const char* sMsg, int nLen, diff --git a/ortools/math_opt/solvers/xpress_solver.cc b/ortools/math_opt/solvers/xpress_solver.cc index 78ec2ee6374..fe214db794f 100644 --- a/ortools/math_opt/solvers/xpress_solver.cc +++ b/ortools/math_opt/solvers/xpress_solver.cc @@ -319,19 +319,18 @@ SolutionStatusProto XpressSolver::getLpSolutionStatus() const { absl::StatusOr> XpressSolver::GetConvexPrimalSolutionIfAvailable( const ModelSolveParametersProto& model_parameters) const { + if (!isFeasible()) { + return SolutionAndClaim{ + .solution = std::nullopt, .feasible_solution_exists = false}; + } PrimalSolutionProto primal_solution; primal_solution.set_feasibility_status(getLpSolutionStatus()); - if (isFeasible()) { - ASSIGN_OR_RETURN(const double sol_val, - xpress_->GetDoubleAttr(XPRS_LPOBJVAL)); - primal_solution.set_objective_value(sol_val); - XpressVectorToSparseDoubleVector(xpress_->GetPrimalValues().value(), - variables_map_, - *primal_solution.mutable_variable_values(), - model_parameters.variable_values_filter()); - } else { - // TODO - } + ASSIGN_OR_RETURN(const double sol_val, xpress_->GetDoubleAttr(XPRS_LPOBJVAL)); + primal_solution.set_objective_value(sol_val); + XpressVectorToSparseDoubleVector(xpress_->GetPrimalValues().value(), + variables_map_, + *primal_solution.mutable_variable_values(), + model_parameters.variable_values_filter()); const bool primal_feasible_solution_exists = (primal_solution.feasibility_status() == SOLUTION_STATUS_FEASIBLE); return SolutionAndClaim{ @@ -342,27 +341,26 @@ XpressSolver::GetConvexPrimalSolutionIfAvailable( absl::StatusOr> XpressSolver::GetConvexDualSolutionIfAvailable( const ModelSolveParametersProto& model_parameters) const { + if (!isFeasible()) { + return SolutionAndClaim{ + .solution = std::nullopt, .feasible_solution_exists = false}; + } DualSolutionProto dual_solution; - const std::vector xprs_constraint_duals = - xpress_->GetConstraintDuals(); + ASSIGN_OR_RETURN(const std::vector xprs_constraint_duals, + xpress_->GetConstraintDuals()); + XpressVectorToSparseDoubleVector(xprs_constraint_duals, linear_constraints_map_, *dual_solution.mutable_dual_values(), model_parameters.dual_values_filter()); - const std::vector xprs_reduced_cost_values = - xpress_->GetReducedCostValues(); + ASSIGN_OR_RETURN(const std::vector xprs_reduced_cost_values, + xpress_->GetReducedCostValues()); XpressVectorToSparseDoubleVector(xprs_reduced_cost_values, variables_map_, *dual_solution.mutable_reduced_costs(), model_parameters.reduced_costs_filter()); - - if (isFeasible()) { - ASSIGN_OR_RETURN(const double sol_val, - xpress_->GetDoubleAttr(XPRS_LPOBJVAL)); - dual_solution.set_objective_value(sol_val); - } else { - // TODO - } + ASSIGN_OR_RETURN(const double sol_val, xpress_->GetDoubleAttr(XPRS_LPOBJVAL)); + dual_solution.set_objective_value(sol_val); dual_solution.set_feasibility_status(getLpSolutionStatus()); bool dual_feasible_solution_exists = (dual_solution.feasibility_status() == SOLUTION_STATUS_FEASIBLE); @@ -397,9 +395,15 @@ inline BasisStatusProto ConvertVariableStatus(const int status) { } absl::StatusOr> XpressSolver::GetBasisIfAvailable() { - // Variable basis BasisProto basis; - auto xprs_variable_basis_status = xpress_->GetVariableBasis(); + + std::vector xprs_variable_basis_status; + std::vector xprs_constraint_basis_status; + RETURN_IF_ERROR(xpress_->GetBasis(xprs_constraint_basis_status, + xprs_variable_basis_status)) + << "Unable to retrieve basis from XPRESS"; + + // Variable basis for (auto [variable_id, xprs_variable_index] : variables_map_) { basis.mutable_variable_status()->add_ids(variable_id); const BasisStatusProto variable_status = @@ -411,12 +415,20 @@ absl::StatusOr> XpressSolver::GetBasisIfAvailable() { } basis.mutable_variable_status()->add_values(variable_status); } + // Constraint basis - // TODO : implement this, mocked for now (else Basis validation will fail) for (auto [constraint_id, xprs_ct_index] : linear_constraints_map_) { basis.mutable_constraint_status()->add_ids(constraint_id); - basis.mutable_constraint_status()->add_values(BASIS_STATUS_BASIC); + const BasisStatusProto status = ConvertVariableStatus( + xprs_constraint_basis_status[xprs_ct_index.constraint_index]); + if (status == BASIS_STATUS_UNSPECIFIED) { + return absl::InternalError(absl::StrCat( + "Invalid Xpress constraint basis status: ", + xprs_constraint_basis_status[xprs_ct_index.constraint_index])); + } + basis.mutable_constraint_status()->add_values(status); } + // Dual basis basis.set_basic_dual_feasibility(SOLUTION_STATUS_UNDETERMINED); if (xpress_status_ == XPRS_LP_OPTIMAL) { @@ -424,6 +436,7 @@ absl::StatusOr> XpressSolver::GetBasisIfAvailable() { } else if (xpress_status_ == XPRS_LP_UNBOUNDED) { basis.set_basic_dual_feasibility(SOLUTION_STATUS_INFEASIBLE); } + return basis; } @@ -454,7 +467,6 @@ void XpressSolver::XpressVectorToSparseDoubleVector( absl::StatusOr XpressSolver::ConvertTerminationReason( SolutionClaims solution_claims, double best_primal_bound, double best_dual_bound) const { - // TODO : improve this if (!is_mip_) { switch (xpress_status_) { case XPRS_LP_UNSTARTED: @@ -476,12 +488,8 @@ absl::StatusOr XpressSolver::ConvertTerminationReason( is_maximize_, LIMIT_UNSPECIFIED, best_primal_bound, best_dual_bound, "Solve did not finish (XPRS_LP_UNFINISHED)"); case XPRS_LP_UNBOUNDED: - if (solution_claims.primal_feasible_solution_exists) { - return UnboundedTerminationProto(is_maximize_); - } - return InfeasibleOrUnboundedTerminationProto( - is_maximize_, FEASIBILITY_STATUS_INFEASIBLE, - "Xpress status XPRS_LP_UNBOUNDED"); + return UnboundedTerminationProto(is_maximize_, + "Xpress status XPRS_LP_UNBOUNDED"); case XPRS_LP_CUTOFF_IN_DUAL: return CutoffTerminationProto( is_maximize_, "Cutoff in dual (XPRS_LP_CUTOFF_IN_DUAL)"); diff --git a/ortools/math_opt/solvers/xpress_solver_test.cc b/ortools/math_opt/solvers/xpress_solver_test.cc index 107bd9e940a..2167bac348a 100644 --- a/ortools/math_opt/solvers/xpress_solver_test.cc +++ b/ortools/math_opt/solvers/xpress_solver_test.cc @@ -32,11 +32,45 @@ using namespace operations_research::math_opt; SimpleLpTestParameters XpressDefaults() { return SimpleLpTestParameters( SolverType::kXpress, SolveParameters(), /*supports_duals=*/true, - /*supports_basis=*/true, + /*supports_basis=*/false, /*ensures_primal_ray=*/false, /*ensures_dual_ray=*/false, /*disallows_infeasible_or_unbounded=*/true); } +TEST(XpressSolverTest, SimpleLpTwoVar) { + double p = 0.0; + auto model = std::make_unique(); + const Variable x_1 = model->AddVariable(0.0, 1.0, false, "x_1"); + const Variable x_2 = model->AddVariable(0.0, 1.0, false, "x_2"); + model->Maximize(2 * x_1 + x_2); + const LinearConstraint y = model->AddLinearConstraint(x_1 + x_2 <= 1.5, "y"); + SolveResult result{Termination::Optimal(2.5)}; + result.solutions.push_back(Solution{ + .primal_solution = + PrimalSolution{.variable_values = {{x_1, 1.0}, {x_2, 0.5}}, + .objective_value = 2.5, + .feasibility_status = SolutionStatus::kFeasible}, + .dual_solution = + DualSolution{.dual_values = {{y, 1.0}}, + .reduced_costs = {{x_1, 1.0}, {x_2, 0.0}}, + .objective_value = 2.5, + .feasibility_status = SolutionStatus::kFeasible}}); + + // Set parameters, e.g. turn on logging. + SolveArguments args; + args.parameters.enable_output = true; + + // Solve and ensure an optimal solution was found with no errors. + const absl::StatusOr resultX = + Solve(*model, SolverType::kXpress, args); + for (auto vb : resultX->solutions[0].basis->variable_status) { + std::cout << vb.first << " = " << vb.second << "\n"; + } + for (auto vb : resultX->solutions[0].basis->constraint_status) { + std::cout << vb.first << " = " << vb.second << "\n"; + } +} + INSTANTIATE_TEST_SUITE_P(XpressSolverLpTest, SimpleLpTest, testing::Values(XpressDefaults())); From d825a91749ed80538be68ce67829d020b04ff489 Mon Sep 17 00:00:00 2001 From: Peter Mitri Date: Wed, 4 Dec 2024 11:40:00 +0100 Subject: [PATCH 07/14] Fix constraint basis. Clean-up test class. --- ortools/math_opt/solvers/xpress/g_xpress.cc | 4 +-- ortools/math_opt/solvers/xpress_solver.cc | 27 +++++++++++--- .../math_opt/solvers/xpress_solver_test.cc | 36 +------------------ 3 files changed, 26 insertions(+), 41 deletions(-) diff --git a/ortools/math_opt/solvers/xpress/g_xpress.cc b/ortools/math_opt/solvers/xpress/g_xpress.cc index 838a69eb878..0ac9facd34d 100644 --- a/ortools/math_opt/solvers/xpress/g_xpress.cc +++ b/ortools/math_opt/solvers/xpress/g_xpress.cc @@ -144,7 +144,7 @@ absl::Status Xpress::ChgCoeffs(absl::Span rowind, } absl::StatusOr Xpress::LpOptimizeAndGetStatus() { - RETURN_IF_ERROR(ToStatus(XPRSlpoptimize(xpress_model_, NULL))) + RETURN_IF_ERROR(ToStatus(XPRSlpoptimize(xpress_model_, nullptr))) << "XPRESS LP solve failed"; int xpress_status; RETURN_IF_ERROR( @@ -157,7 +157,7 @@ absl::Status Xpress::PostSolve() { } absl::StatusOr Xpress::MipOptimizeAndGetStatus() { - RETURN_IF_ERROR(ToStatus(XPRSmipoptimize(xpress_model_, NULL))) + RETURN_IF_ERROR(ToStatus(XPRSmipoptimize(xpress_model_, nullptr))) << "XPRESS MIP solve failed"; int xpress_status; RETURN_IF_ERROR( diff --git a/ortools/math_opt/solvers/xpress_solver.cc b/ortools/math_opt/solvers/xpress_solver.cc index fe214db794f..2f7f7b18ba4 100644 --- a/ortools/math_opt/solvers/xpress_solver.cc +++ b/ortools/math_opt/solvers/xpress_solver.cc @@ -394,6 +394,23 @@ inline BasisStatusProto ConvertVariableStatus(const int status) { } } +inline BasisStatusProto ConvertConstraintStatus(const int status) { + // XPRESS row basis status is that of the slack variable + // For example, if the slack variable is at LB, the constraint is at UB + switch (status) { + case XPRS_BASIC: + return BASIS_STATUS_BASIC; + case XPRS_AT_LOWER: + return BASIS_STATUS_AT_UPPER_BOUND; + case XPRS_AT_UPPER: + return BASIS_STATUS_AT_LOWER_BOUND; + case XPRS_FREE_SUPER: + return BASIS_STATUS_FREE; + default: + return BASIS_STATUS_UNSPECIFIED; + } +} + absl::StatusOr> XpressSolver::GetBasisIfAvailable() { BasisProto basis; @@ -419,7 +436,7 @@ absl::StatusOr> XpressSolver::GetBasisIfAvailable() { // Constraint basis for (auto [constraint_id, xprs_ct_index] : linear_constraints_map_) { basis.mutable_constraint_status()->add_ids(constraint_id); - const BasisStatusProto status = ConvertVariableStatus( + const BasisStatusProto status = ConvertConstraintStatus( xprs_constraint_basis_status[xprs_ct_index.constraint_index]); if (status == BASIS_STATUS_UNSPECIFIED) { return absl::InternalError(absl::StrCat( @@ -513,9 +530,11 @@ absl::StatusOr XpressSolver::ConvertTerminationReason( absl::StatusOr XpressSolver::Update( const ModelUpdateProto& model_update) { - // TODO: implement this - return absl::UnimplementedError( - "XpressSolver::Update is not implemented yet"); + if (!UpdateIsSupported(model_update, kXpressSupportedStructures)) { + return false; + } + // TODO: implement Update + return false; } absl::StatusOr diff --git a/ortools/math_opt/solvers/xpress_solver_test.cc b/ortools/math_opt/solvers/xpress_solver_test.cc index 2167bac348a..107bd9e940a 100644 --- a/ortools/math_opt/solvers/xpress_solver_test.cc +++ b/ortools/math_opt/solvers/xpress_solver_test.cc @@ -32,45 +32,11 @@ using namespace operations_research::math_opt; SimpleLpTestParameters XpressDefaults() { return SimpleLpTestParameters( SolverType::kXpress, SolveParameters(), /*supports_duals=*/true, - /*supports_basis=*/false, + /*supports_basis=*/true, /*ensures_primal_ray=*/false, /*ensures_dual_ray=*/false, /*disallows_infeasible_or_unbounded=*/true); } -TEST(XpressSolverTest, SimpleLpTwoVar) { - double p = 0.0; - auto model = std::make_unique(); - const Variable x_1 = model->AddVariable(0.0, 1.0, false, "x_1"); - const Variable x_2 = model->AddVariable(0.0, 1.0, false, "x_2"); - model->Maximize(2 * x_1 + x_2); - const LinearConstraint y = model->AddLinearConstraint(x_1 + x_2 <= 1.5, "y"); - SolveResult result{Termination::Optimal(2.5)}; - result.solutions.push_back(Solution{ - .primal_solution = - PrimalSolution{.variable_values = {{x_1, 1.0}, {x_2, 0.5}}, - .objective_value = 2.5, - .feasibility_status = SolutionStatus::kFeasible}, - .dual_solution = - DualSolution{.dual_values = {{y, 1.0}}, - .reduced_costs = {{x_1, 1.0}, {x_2, 0.0}}, - .objective_value = 2.5, - .feasibility_status = SolutionStatus::kFeasible}}); - - // Set parameters, e.g. turn on logging. - SolveArguments args; - args.parameters.enable_output = true; - - // Solve and ensure an optimal solution was found with no errors. - const absl::StatusOr resultX = - Solve(*model, SolverType::kXpress, args); - for (auto vb : resultX->solutions[0].basis->variable_status) { - std::cout << vb.first << " = " << vb.second << "\n"; - } - for (auto vb : resultX->solutions[0].basis->constraint_status) { - std::cout << vb.first << " = " << vb.second << "\n"; - } -} - INSTANTIATE_TEST_SUITE_P(XpressSolverLpTest, SimpleLpTest, testing::Values(XpressDefaults())); From d0502eaac069898f729f3cb741dc2b1673f66c1c Mon Sep 17 00:00:00 2001 From: Peter Mitri Date: Thu, 5 Dec 2024 11:41:55 +0100 Subject: [PATCH 08/14] Activate some test suites and improve XpressSolver error handling --- ortools/math_opt/solvers/CMakeLists.txt | 8 +- ortools/math_opt/solvers/xpress/g_xpress.cc | 45 ++++++--- ortools/math_opt/solvers/xpress/g_xpress.h | 7 +- ortools/math_opt/solvers/xpress_solver.cc | 78 ++++++++++++--- ortools/math_opt/solvers/xpress_solver.h | 3 + .../math_opt/solvers/xpress_solver_test.cc | 95 ++++++++++++++----- 6 files changed, 181 insertions(+), 55 deletions(-) diff --git a/ortools/math_opt/solvers/CMakeLists.txt b/ortools/math_opt/solvers/CMakeLists.txt index f6af550c765..64cebe10015 100644 --- a/ortools/math_opt/solvers/CMakeLists.txt +++ b/ortools/math_opt/solvers/CMakeLists.txt @@ -251,17 +251,21 @@ if(USE_XPRESS) absl::status ortools::math_opt_matchers "$" + "$" "$" "$" "$" "$" - #"$" - #"$" + "$" + "$" "$" "$" "$" "$" "$" "$" + "$" + "$" + "$" ) endif() diff --git a/ortools/math_opt/solvers/xpress/g_xpress.cc b/ortools/math_opt/solvers/xpress/g_xpress.cc index 0ac9facd34d..2eb335169fe 100644 --- a/ortools/math_opt/solvers/xpress/g_xpress.cc +++ b/ortools/math_opt/solvers/xpress/g_xpress.cc @@ -187,52 +187,71 @@ absl::StatusOr Xpress::GetDoubleAttr(int attribute) const { } absl::StatusOr> Xpress::GetPrimalValues() const { - int nCols = GetNumberOfColumns(); - XPRSgetintattrib(xpress_model_, XPRS_COLS, &nCols); - std::vector values(nCols); + int nVars = GetNumberOfVariables(); + XPRSgetintattrib(xpress_model_, XPRS_COLS, &nVars); + std::vector values(nVars); RETURN_IF_ERROR(ToStatus( XPRSgetlpsol(xpress_model_, values.data(), nullptr, nullptr, nullptr))) << "Error getting Xpress LP solution"; return values; } -int Xpress::GetNumberOfRows() const { +int Xpress::GetNumberOfConstraints() const { int n; XPRSgetintattrib(xpress_model_, XPRS_ROWS, &n); return n; } -int Xpress::GetNumberOfColumns() const { +int Xpress::GetNumberOfVariables() const { int n; XPRSgetintattrib(xpress_model_, XPRS_COLS, &n); return n; } absl::StatusOr> Xpress::GetConstraintDuals() const { - int nRows = GetNumberOfRows(); - double values[nRows]; + int nCons = GetNumberOfConstraints(); + double values[nCons]; RETURN_IF_ERROR( ToStatus(XPRSgetlpsol(xpress_model_, nullptr, nullptr, values, nullptr))) << "Failed to retrieve LP solution from XPRESS"; - std::vector result(values, values + nRows); + std::vector result(values, values + nCons); return result; } absl::StatusOr> Xpress::GetReducedCostValues() const { - int nCols = GetNumberOfColumns(); - double values[nCols]; + int nVars = GetNumberOfVariables(); + double values[nVars]; RETURN_IF_ERROR( ToStatus(XPRSgetlpsol(xpress_model_, nullptr, nullptr, nullptr, values))) << "Failed to retrieve LP solution from XPRESS"; - std::vector result(values, values + nCols); + std::vector result(values, values + nVars); return result; } absl::Status Xpress::GetBasis(std::vector& rowBasis, std::vector& colBasis) const { - rowBasis.resize(GetNumberOfRows()); - colBasis.resize(GetNumberOfColumns()); + rowBasis.resize(GetNumberOfConstraints()); + colBasis.resize(GetNumberOfVariables()); return ToStatus( XPRSgetbasis(xpress_model_, rowBasis.data(), colBasis.data())); } +absl::StatusOr> Xpress::GetVarLb() const { + int nVars = GetNumberOfVariables(); + std::vector bounds; + bounds.reserve(nVars); + RETURN_IF_ERROR( + ToStatus(XPRSgetlb(xpress_model_, bounds.data(), 0, nVars - 1))) + << "Failed to retrieve variable LB from XPRESS"; + return bounds; +} +absl::StatusOr> Xpress::GetVarUb() const { + int nVars = GetNumberOfVariables(); + std::vector bounds; + bounds.reserve(nVars); + RETURN_IF_ERROR( + ToStatus(XPRSgetub(xpress_model_, bounds.data(), 0, nVars - 1))) + << "Failed to retrieve variable UB from XPRESS"; + return bounds; +} + } // namespace operations_research::math_opt diff --git a/ortools/math_opt/solvers/xpress/g_xpress.h b/ortools/math_opt/solvers/xpress/g_xpress.h index 05419536256..d7fc9c1aa7c 100644 --- a/ortools/math_opt/solvers/xpress/g_xpress.h +++ b/ortools/math_opt/solvers/xpress/g_xpress.h @@ -89,8 +89,11 @@ class Xpress { const char* sMsg, int nLen, int nMsgLvl); - int GetNumberOfRows() const; - int GetNumberOfColumns() const; + int GetNumberOfConstraints() const; + int GetNumberOfVariables() const; + + absl::StatusOr> GetVarLb() const; + absl::StatusOr> GetVarUb() const; private: XPRSprob xpress_model_; diff --git a/ortools/math_opt/solvers/xpress_solver.cc b/ortools/math_opt/solvers/xpress_solver.cc index 2f7f7b18ba4..1138d0e942f 100644 --- a/ortools/math_opt/solvers/xpress_solver.cc +++ b/ortools/math_opt/solvers/xpress_solver.cc @@ -13,14 +13,33 @@ #include "ortools/math_opt/solvers/xpress_solver.h" +#include "absl/strings/str_join.h" #include "ortools/base/map_util.h" #include "ortools/base/protoutil.h" #include "ortools/base/status_macros.h" #include "ortools/math_opt/core/math_opt_proto_utils.h" #include "ortools/math_opt/cpp/solve_result.h" +#include "ortools/math_opt/validators/callback_validator.h" #include "ortools/xpress/environment.h" -namespace operations_research::math_opt { +namespace operations_research { +namespace math_opt { +namespace { +// TODO: set parameters +// TODO: add all unsupported parameter combinations here +absl::Status SetParameters(const SolveParametersProto& parameters) { + std::vector warnings; + if (parameters.has_threads() && parameters.threads() > 1) { + warnings.push_back(absl::StrCat( + "XpressSolver only supports parameters.threads = 1; value ", + parameters.threads(), " is not supported")); + } + if (!warnings.empty()) { + return absl::InvalidArgumentError(absl::StrJoin(warnings, "; ")); + } + return absl::OkStatus(); +} +} // namespace constexpr SupportedProblemStructures kXpressSupportedStructures = { .integer_variables = SupportType::kNotImplemented, @@ -66,7 +85,7 @@ absl::Status XpressSolver::AddNewVariables( const VariablesProto& new_variables) { const int num_new_variables = new_variables.lower_bounds().size(); std::vector variable_type(num_new_variables); - int n_variables = xpress_->GetNumberOfColumns(); + int n_variables = xpress_->GetNumberOfVariables(); for (int j = 0; j < num_new_variables; ++j) { const VarId id = new_variables.ids(j); InsertOrDie(&variables_map_, id, j + n_variables); @@ -101,7 +120,7 @@ absl::Status XpressSolver::AddNewLinearConstraints( constraint_rhs.reserve(num_new_constraints); constraint_sense.reserve(num_new_constraints); new_slacks.reserve(num_new_constraints); - int n_constraints = xpress_->GetNumberOfRows(); + int n_constraints = xpress_->GetNumberOfConstraints(); for (int i = 0; i < num_new_constraints; ++i) { const int64_t id = constraints.ids(i); LinearConstraintData& constraint_data = @@ -126,9 +145,6 @@ absl::Status XpressSolver::AddNewLinearConstraints( sense = XPRS_EQUAL; rhs = lb; } else { - if (ub < lb) { - return absl::InvalidArgumentError("Lower bound > Upper bound"); - } sense = XPRS_RANGE; rhs = ub; rng = ub - lb; @@ -149,9 +165,6 @@ absl::Status XpressSolver::AddSingleObjective(const ObjectiveProto& objective) { "Quadratic objectives are not yet implemented in XPRESS solver " "interface."); }*/ - if (objective.linear_coefficients().ids_size() == 0) { - return absl::OkStatus(); - } std::vector index; index.reserve(objective.linear_coefficients().ids_size()); for (const int64_t id : objective.linear_coefficients().ids()) { @@ -194,6 +207,18 @@ absl::StatusOr XpressSolver::Solve( // TODO: set lazy constraints // TODO: add interrupter using xpress_->Terminate(); + RETURN_IF_ERROR(CheckRegisteredCallbackEvents(callback_registration, + /*supported_events=*/{})); + + RETURN_IF_ERROR(SetParameters(parameters)); + + // Check that bounds are not inverted just before solve + // XPRESS returns "infeasible" when bounds are inverted + { + ASSIGN_OR_RETURN(const InvertedBounds inv_bounds, ListInvertedBounds()); + RETURN_IF_ERROR(inv_bounds.ToStatus()); + } + RETURN_IF_ERROR(CallXpressSolve(parameters.enable_output())) << "Error during XPRESS solve"; @@ -325,7 +350,7 @@ XpressSolver::GetConvexPrimalSolutionIfAvailable( } PrimalSolutionProto primal_solution; primal_solution.set_feasibility_status(getLpSolutionStatus()); - ASSIGN_OR_RETURN(const double sol_val, xpress_->GetDoubleAttr(XPRS_LPOBJVAL)); + ASSIGN_OR_RETURN(const double sol_val, GetBestPrimalBound()); primal_solution.set_objective_value(sol_val); XpressVectorToSparseDoubleVector(xpress_->GetPrimalValues().value(), variables_map_, @@ -359,13 +384,12 @@ XpressSolver::GetConvexDualSolutionIfAvailable( XpressVectorToSparseDoubleVector(xprs_reduced_cost_values, variables_map_, *dual_solution.mutable_reduced_costs(), model_parameters.reduced_costs_filter()); - ASSIGN_OR_RETURN(const double sol_val, xpress_->GetDoubleAttr(XPRS_LPOBJVAL)); + ASSIGN_OR_RETURN(const double sol_val, GetBestPrimalBound()); dual_solution.set_objective_value(sol_val); dual_solution.set_feasibility_status(getLpSolutionStatus()); bool dual_feasible_solution_exists = (dual_solution.feasibility_status() == SOLUTION_STATUS_FEASIBLE); - ASSIGN_OR_RETURN(const double best_dual_bound, - xpress_->GetDoubleAttr(XPRS_LPOBJVAL)); + ASSIGN_OR_RETURN(const double best_dual_bound, GetBestDualBound()); // const double best_dual_bound = is_maximize_ ? kMinusInf : kPlusInf; if (dual_feasible_solution_exists || std::isfinite(best_dual_bound)) { dual_feasible_solution_exists = true; @@ -546,6 +570,30 @@ XpressSolver::ComputeInfeasibleSubsystem(const SolveParametersProto& parameters, "XpressSolver::ComputeInfeasibleSubsystem is not implemented yet"); } -MATH_OPT_REGISTER_SOLVER(SOLVER_TYPE_XPRESS, XpressSolver::New) +absl::StatusOr XpressSolver::ListInvertedBounds() const { + InvertedBounds inverted_bounds; + { + ASSIGN_OR_RETURN(const std::vector var_lbs, xpress_->GetVarLb()); + ASSIGN_OR_RETURN(const std::vector var_ubs, xpress_->GetVarUb()); + for (const auto& [id, index] : variables_map_) { + if (var_lbs[index] > var_ubs[index]) { + inverted_bounds.variables.push_back(id); + } + } + } + // TODO: we better fetch constraint lb & ub from xpress_ to ensure sync + for (const auto& [id, cstr_data] : linear_constraints_map_) { + if (cstr_data.lower_bound > cstr_data.upper_bound) { + inverted_bounds.linear_constraints.push_back(id); + } + } + // Above code have inserted ids in non-stable order. + std::sort(inverted_bounds.variables.begin(), inverted_bounds.variables.end()); + std::sort(inverted_bounds.linear_constraints.begin(), + inverted_bounds.linear_constraints.end()); + return inverted_bounds; +} -} // namespace operations_research::math_opt +MATH_OPT_REGISTER_SOLVER(SOLVER_TYPE_XPRESS, XpressSolver::New) +} // namespace math_opt +} // namespace operations_research diff --git a/ortools/math_opt/solvers/xpress_solver.h b/ortools/math_opt/solvers/xpress_solver.h index 5474d829da5..62899ecad3e 100644 --- a/ortools/math_opt/solvers/xpress_solver.h +++ b/ortools/math_opt/solvers/xpress_solver.h @@ -26,6 +26,7 @@ #include "absl/types/span.h" #include "ortools/base/linked_hash_map.h" #include "ortools/math_opt/callback.pb.h" +#include "ortools/math_opt/core/inverted_bounds.h" #include "ortools/math_opt/core/solver_interface.h" #include "ortools/math_opt/infeasible_subsystem.pb.h" #include "ortools/math_opt/model.pb.h" @@ -195,6 +196,8 @@ class XpressSolver : public SolverInterface { SolutionStatusProto getLpSolutionStatus() const; + absl::StatusOr ListInvertedBounds() const; + // Fields to track the number of Xpress variables and constraints. These // quantities are updated immediately after adding or removing to the model, // so it is correct even if XPRESS C API has not yet been called. diff --git a/ortools/math_opt/solvers/xpress_solver_test.cc b/ortools/math_opt/solvers/xpress_solver_test.cc index 107bd9e940a..7c2150089e6 100644 --- a/ortools/math_opt/solvers/xpress_solver_test.cc +++ b/ortools/math_opt/solvers/xpress_solver_test.cc @@ -11,36 +11,85 @@ // See the License for the specific language governing permissions and // limitations under the License. -#include +#include +#include +#include +#include -#include - -#include "absl/log/check.h" -#include "absl/status/statusor.h" #include "gtest/gtest.h" -#include "ortools/base/gmock.h" -#include "ortools/base/init_google.h" -#include "ortools/math_opt/core/solver.h" -#include "ortools/math_opt/cpp/matchers.h" #include "ortools/math_opt/cpp/math_opt.h" +#include "ortools/math_opt/solver_tests/callback_tests.h" +#include "ortools/math_opt/solver_tests/generic_tests.h" +#include "ortools/math_opt/solver_tests/infeasible_subsystem_tests.h" +#include "ortools/math_opt/solver_tests/invalid_input_tests.h" +#include "ortools/math_opt/solver_tests/logical_constraint_tests.h" +#include "ortools/math_opt/solver_tests/lp_incomplete_solve_tests.h" +#include "ortools/math_opt/solver_tests/lp_initial_basis_tests.h" #include "ortools/math_opt/solver_tests/lp_model_solve_parameters_tests.h" +#include "ortools/math_opt/solver_tests/lp_parameter_tests.h" #include "ortools/math_opt/solver_tests/lp_tests.h" -#include "ortools/port/proto_utils.h" +#include "ortools/math_opt/solver_tests/multi_objective_tests.h" +#include "ortools/math_opt/solver_tests/qc_tests.h" +#include "ortools/math_opt/solver_tests/qp_tests.h" +#include "ortools/math_opt/solver_tests/second_order_cone_tests.h" +#include "ortools/math_opt/solver_tests/status_tests.h" -using namespace operations_research::math_opt; +namespace operations_research { +namespace math_opt { +namespace { +using testing::ValuesIn; + +// TODO: implement missing LP features +INSTANTIATE_TEST_SUITE_P( + XpressSolverLpTest, SimpleLpTest, + testing::Values(SimpleLpTestParameters( + SolverType::kXpress, SolveParameters(), /*supports_duals=*/true, + /*supports_basis=*/true, + /*ensures_primal_ray=*/false, /*ensures_dual_ray=*/false, + /*disallows_infeasible_or_unbounded=*/true))); + +// TODO: implement message callbacks +INSTANTIATE_TEST_SUITE_P(XpressMessageCallbackTest, MessageCallbackTest, + testing::Values(MessageCallbackTestParams( + SolverType::kXpress, + /*support_message_callback=*/false, + /*support_interrupter=*/false, + /*integer_variables=*/false, ""))); -SimpleLpTestParameters XpressDefaults() { - return SimpleLpTestParameters( - SolverType::kXpress, SolveParameters(), /*supports_duals=*/true, - /*supports_basis=*/true, - /*ensures_primal_ray=*/false, /*ensures_dual_ray=*/false, - /*disallows_infeasible_or_unbounded=*/true); -} -INSTANTIATE_TEST_SUITE_P(XpressSolverLpTest, SimpleLpTest, - testing::Values(XpressDefaults())); +// TODO: implement callbacks +INSTANTIATE_TEST_SUITE_P( + XpressCallbackTest, CallbackTest, + testing::Values(CallbackTestParams(SolverType::kXpress, + /*integer_variables=*/false, + /*add_lazy_constraints=*/false, + /*add_cuts=*/false, + /*supported_events=*/{}, + /*all_solutions=*/std::nullopt, + /*reaches_cut_callback*/ std::nullopt))); -int main(int argc, char** argv) { - testing::InitGoogleTest(&argc, argv); - RUN_ALL_TESTS(); +INSTANTIATE_TEST_SUITE_P(XpressInvalidInputTest, InvalidInputTest, + testing::Values(InvalidInputTestParameters( + SolverType::kXpress, + /*use_integer_variables=*/false))); + +InvalidParameterTestParams InvalidThreadsParameters() { + SolveParameters params; + params.threads = 2; + return InvalidParameterTestParams(SolverType::kXpress, std::move(params), + {"only supports parameters.threads = 1"}); } + +// TODO: add all invalid parameters combinations +INSTANTIATE_TEST_SUITE_P(XpressInvalidParameterTest, InvalidParameterTest, +ValuesIn({InvalidThreadsParameters()})); + +INSTANTIATE_TEST_SUITE_P(XpressGenericTest, GenericTest, + testing::Values(GenericTestParameters( + SolverType::kXpress, /*support_interrupter=*/false, + /*integer_variables=*/false, + /*expected_log=*/"Optimal solution found"))); + +} // namespace +} // namespace math_opt +} // namespace operations_research \ No newline at end of file From dc317c07a84599a60346757437aaa77542089e44 Mon Sep 17 00:00:00 2001 From: Peter Mitri Date: Thu, 5 Dec 2024 11:59:09 +0100 Subject: [PATCH 09/14] Add IIS test --- ortools/math_opt/cpp/parameters.cc | 2 ++ ortools/math_opt/solvers/xpress_solver.cc | 2 +- ortools/math_opt/solvers/xpress_solver_test.cc | 13 +++++++++++-- ortools/service/v1/mathopt/parameters.proto | 6 ++++++ 4 files changed, 20 insertions(+), 3 deletions(-) diff --git a/ortools/math_opt/cpp/parameters.cc b/ortools/math_opt/cpp/parameters.cc index 562367beae8..4469168bb79 100644 --- a/ortools/math_opt/cpp/parameters.cc +++ b/ortools/math_opt/cpp/parameters.cc @@ -85,6 +85,8 @@ std::optional Enum::ToOptString( return "highs"; case SolverType::kSantorini: return "santorini"; + case SolverType::kXpress: + return "xpress"; } return std::nullopt; } diff --git a/ortools/math_opt/solvers/xpress_solver.cc b/ortools/math_opt/solvers/xpress_solver.cc index 1138d0e942f..beb8128d1bc 100644 --- a/ortools/math_opt/solvers/xpress_solver.cc +++ b/ortools/math_opt/solvers/xpress_solver.cc @@ -567,7 +567,7 @@ XpressSolver::ComputeInfeasibleSubsystem(const SolveParametersProto& parameters, const SolveInterrupter* interrupter) { // TODO: implement this return absl::UnimplementedError( - "XpressSolver::ComputeInfeasibleSubsystem is not implemented yet"); + "XpressSolver cannot compute infeasible subsystem yet"); } absl::StatusOr XpressSolver::ListInvertedBounds() const { diff --git a/ortools/math_opt/solvers/xpress_solver_test.cc b/ortools/math_opt/solvers/xpress_solver_test.cc index 7c2150089e6..f038e3215cf 100644 --- a/ortools/math_opt/solvers/xpress_solver_test.cc +++ b/ortools/math_opt/solvers/xpress_solver_test.cc @@ -56,7 +56,6 @@ INSTANTIATE_TEST_SUITE_P(XpressMessageCallbackTest, MessageCallbackTest, /*support_interrupter=*/false, /*integer_variables=*/false, ""))); - // TODO: implement callbacks INSTANTIATE_TEST_SUITE_P( XpressCallbackTest, CallbackTest, @@ -82,7 +81,7 @@ InvalidParameterTestParams InvalidThreadsParameters() { // TODO: add all invalid parameters combinations INSTANTIATE_TEST_SUITE_P(XpressInvalidParameterTest, InvalidParameterTest, -ValuesIn({InvalidThreadsParameters()})); + ValuesIn({InvalidThreadsParameters()})); INSTANTIATE_TEST_SUITE_P(XpressGenericTest, GenericTest, testing::Values(GenericTestParameters( @@ -90,6 +89,16 @@ INSTANTIATE_TEST_SUITE_P(XpressGenericTest, GenericTest, /*integer_variables=*/false, /*expected_log=*/"Optimal solution found"))); +// TODO: When XPRESS callbacks are supported, enable this test. +GTEST_ALLOW_UNINSTANTIATED_PARAMETERIZED_TEST(TimeLimitTest); + +// TODO: implement IIS support +INSTANTIATE_TEST_SUITE_P(XpressInfeasibleSubsystemTest, InfeasibleSubsystemTest, + testing::Values(InfeasibleSubsystemTestParameters( + {.solver_type = SolverType::kXpress, + .support_menu = { + .supports_infeasible_subsystems = false}}))); + } // namespace } // namespace math_opt } // namespace operations_research \ No newline at end of file diff --git a/ortools/service/v1/mathopt/parameters.proto b/ortools/service/v1/mathopt/parameters.proto index 5dc0e4e524c..99a57ee5318 100644 --- a/ortools/service/v1/mathopt/parameters.proto +++ b/ortools/service/v1/mathopt/parameters.proto @@ -99,6 +99,12 @@ enum SolverTypeProto { // Slow/not recommended for production. Not an LP solver (no dual information // returned). SOLVER_TYPE_SANTORINI = 11; + + // Fico XPRESS solver (third party). + // + // Supports LP, MIP, and nonconvex integer quadratic problems. + // A fast option, but has special licensing. + SOLVER_TYPE_XPRESS = 12; } // Selects an algorithm for solving linear programs. From 62e6b26fb61b1caa61a58a210c0ca539f9a7ab24 Mon Sep 17 00:00:00 2001 From: Peter Mitri Date: Thu, 5 Dec 2024 13:12:03 +0100 Subject: [PATCH 10/14] Enable LpParameterTestParams and add LP parameters support (algorithm & iter limits) --- .../math_opt/solver_tests/base_solver_test.cc | 6 + ortools/math_opt/solvers/xpress/g_xpress.cc | 4 +- ortools/math_opt/solvers/xpress/g_xpress.h | 2 +- ortools/math_opt/solvers/xpress_solver.cc | 111 ++++++++++++++++-- ortools/math_opt/solvers/xpress_solver.h | 2 +- .../math_opt/solvers/xpress_solver_test.cc | 34 +++++- ortools/xpress/environment.h | 1 + 7 files changed, 142 insertions(+), 18 deletions(-) diff --git a/ortools/math_opt/solver_tests/base_solver_test.cc b/ortools/math_opt/solver_tests/base_solver_test.cc index f88a9084016..1efb52b8c4d 100644 --- a/ortools/math_opt/solver_tests/base_solver_test.cc +++ b/ortools/math_opt/solver_tests/base_solver_test.cc @@ -50,6 +50,9 @@ bool ActivatePrimalRay(const SolverType solver_type, SolveParameters& params) { return false; case SolverType::kHighs: return false; + case SolverType::kXpress: + // TODO: support XPRESS + return false; default: LOG(FATAL) << "Solver " << solver_type @@ -82,6 +85,9 @@ bool ActivateDualRay(const SolverType solver_type, SolveParameters& params) { return false; case SolverType::kHighs: return false; + case SolverType::kXpress: + // TODO: support XPRESS + return false; default: LOG(FATAL) << "Solver " << solver_type diff --git a/ortools/math_opt/solvers/xpress/g_xpress.cc b/ortools/math_opt/solvers/xpress/g_xpress.cc index 2eb335169fe..00d7baf0390 100644 --- a/ortools/math_opt/solvers/xpress/g_xpress.cc +++ b/ortools/math_opt/solvers/xpress/g_xpress.cc @@ -143,8 +143,8 @@ absl::Status Xpress::ChgCoeffs(absl::Span rowind, XPRSchgmcoef(xpress_model_, n_coefs, c_rowind, c_colind, c_values)); } -absl::StatusOr Xpress::LpOptimizeAndGetStatus() { - RETURN_IF_ERROR(ToStatus(XPRSlpoptimize(xpress_model_, nullptr))) +absl::StatusOr Xpress::LpOptimizeAndGetStatus(std::string flags) { + RETURN_IF_ERROR(ToStatus(XPRSlpoptimize(xpress_model_, flags.c_str()))) << "XPRESS LP solve failed"; int xpress_status; RETURN_IF_ERROR( diff --git a/ortools/math_opt/solvers/xpress/g_xpress.h b/ortools/math_opt/solvers/xpress/g_xpress.h index d7fc9c1aa7c..83b5cddb6c6 100644 --- a/ortools/math_opt/solvers/xpress/g_xpress.h +++ b/ortools/math_opt/solvers/xpress/g_xpress.h @@ -74,7 +74,7 @@ class Xpress { absl::Status ChgCoeffs(absl::Span cind, absl::Span vind, absl::Span val); - absl::StatusOr LpOptimizeAndGetStatus(); + absl::StatusOr LpOptimizeAndGetStatus(std::string flags); absl::StatusOr MipOptimizeAndGetStatus(); absl::Status PostSolve(); diff --git a/ortools/math_opt/solvers/xpress_solver.cc b/ortools/math_opt/solvers/xpress_solver.cc index beb8128d1bc..6f11b0633db 100644 --- a/ortools/math_opt/solvers/xpress_solver.cc +++ b/ortools/math_opt/solvers/xpress_solver.cc @@ -20,20 +20,37 @@ #include "ortools/math_opt/core/math_opt_proto_utils.h" #include "ortools/math_opt/cpp/solve_result.h" #include "ortools/math_opt/validators/callback_validator.h" +#include "ortools/port/proto_utils.h" #include "ortools/xpress/environment.h" namespace operations_research { namespace math_opt { namespace { -// TODO: set parameters // TODO: add all unsupported parameter combinations here -absl::Status SetParameters(const SolveParametersProto& parameters) { +absl::Status CheckParameters(const SolveParametersProto& parameters) { std::vector warnings; if (parameters.has_threads() && parameters.threads() > 1) { warnings.push_back(absl::StrCat( "XpressSolver only supports parameters.threads = 1; value ", parameters.threads(), " is not supported")); } + if (parameters.lp_algorithm() != LP_ALGORITHM_UNSPECIFIED && + parameters.lp_algorithm() != LP_ALGORITHM_PRIMAL_SIMPLEX && + parameters.lp_algorithm() != LP_ALGORITHM_DUAL_SIMPLEX && + parameters.lp_algorithm() != LP_ALGORITHM_BARRIER) { + warnings.emplace_back(absl::StrCat( + "XpressSolver does not support the 'lp_algorithm' parameter value: ", + ProtoEnumToString(parameters.lp_algorithm()))); + } + if (parameters.has_objective_limit()) { + warnings.emplace_back("XpressSolver does not support objective_limit yet"); + } + if (parameters.has_best_bound_limit()) { + warnings.emplace_back("XpressSolver does not support best_bound_limit yet"); + } + if (parameters.has_cutoff_limit()) { + warnings.emplace_back("XpressSolver does not support cutoff_limit yet"); + } if (!warnings.empty()) { return absl::InvalidArgumentError(absl::StrJoin(warnings, "; ")); } @@ -210,7 +227,7 @@ absl::StatusOr XpressSolver::Solve( RETURN_IF_ERROR(CheckRegisteredCallbackEvents(callback_registration, /*supported_events=*/{})); - RETURN_IF_ERROR(SetParameters(parameters)); + RETURN_IF_ERROR(CheckParameters(parameters)); // Check that bounds are not inverted just before solve // XPRESS returns "infeasible" when bounds are inverted @@ -219,8 +236,7 @@ absl::StatusOr XpressSolver::Solve( RETURN_IF_ERROR(inv_bounds.ToStatus()); } - RETURN_IF_ERROR(CallXpressSolve(parameters.enable_output())) - << "Error during XPRESS solve"; + RETURN_IF_ERROR(CallXpressSolve(parameters)) << "Error during XPRESS solve"; ASSIGN_OR_RETURN(SolveResultProto solve_result, ExtractSolveResultProto(start, model_parameters)); @@ -228,9 +244,10 @@ absl::StatusOr XpressSolver::Solve( return solve_result; } -absl::Status XpressSolver::CallXpressSolve(bool enableOutput) { +absl::Status XpressSolver::CallXpressSolve( + const SolveParametersProto& parameters) { // Enable screen output right before solve - if (enableOutput) { + if (parameters.enable_output()) { RETURN_IF_ERROR(xpress_->SetIntAttr(XPRS_OUTPUTLOG, 1)) << "Unable to enable XPRESS logs"; } @@ -238,7 +255,50 @@ absl::Status XpressSolver::CallXpressSolve(bool enableOutput) { if (is_mip_) { ASSIGN_OR_RETURN(xpress_status_, xpress_->MipOptimizeAndGetStatus()); } else { - ASSIGN_OR_RETURN(xpress_status_, xpress_->LpOptimizeAndGetStatus()); + std::string flags; + // TODO: put "p", "d", "b" in environment.h? + // TODO: refactor this, it is ugly + switch (parameters.lp_algorithm()) { + case LP_ALGORITHM_PRIMAL_SIMPLEX: + flags = "p"; + if (parameters.has_iteration_limit()) { + RETURN_IF_ERROR(xpress_->SetIntAttr(XPRS_LPITERLIMIT, + parameters.iteration_limit())); + } + break; + case LP_ALGORITHM_DUAL_SIMPLEX: + flags = "d"; + if (parameters.has_iteration_limit()) { + RETURN_IF_ERROR(xpress_->SetIntAttr(XPRS_LPITERLIMIT, + parameters.iteration_limit())); + } + break; + case LP_ALGORITHM_BARRIER: + flags = "b"; + if (parameters.has_iteration_limit()) { + RETURN_IF_ERROR(xpress_->SetIntAttr(XPRS_BARITERLIMIT, + parameters.iteration_limit())); + // we also have to set simplex iter limits because barrier iterations + // call simplex iterations + // TODO: verify this + RETURN_IF_ERROR(xpress_->SetIntAttr(XPRS_LPITERLIMIT, + parameters.iteration_limit())); + } + break; + default: + // this makes XPRESS use default algorithm (XPRS_DEFAULTALG) + flags = ""; + if (parameters.has_iteration_limit()) { + // TODO: fetch the default alg and set its limit? or default to primal + // simplex + RETURN_IF_ERROR(xpress_->SetIntAttr(XPRS_LPITERLIMIT, + parameters.iteration_limit())); + RETURN_IF_ERROR(xpress_->SetIntAttr(XPRS_BARITERLIMIT, + parameters.iteration_limit())); + } + break; + } + ASSIGN_OR_RETURN(xpress_status_, xpress_->LpOptimizeAndGetStatus(flags)); } // Post-solve if (!(is_mip_ ? (xpress_status_ == XPRS_MIP_OPTIMAL) @@ -246,7 +306,7 @@ absl::Status XpressSolver::CallXpressSolve(bool enableOutput) { RETURN_IF_ERROR(xpress_->PostSolve()) << "Post-solve failed in XPRESS"; } // Disable screen output right after solve - if (enableOutput) { + if (parameters.enable_output()) { RETURN_IF_ERROR(xpress_->SetIntAttr(XPRS_OUTPUTLOG, 0)) << "Unable to disable XPRESS logs"; } @@ -324,17 +384,27 @@ absl::StatusOr XpressSolver::GetLpSolution( } bool XpressSolver::isFeasible() const { - return xpress_status_ == (is_mip_ ? XPRS_MIP_OPTIMAL : XPRS_LP_OPTIMAL); + if (is_mip_) { + return xpress_status_ == XPRS_MIP_OPTIMAL; + } else { + return xpress_status_ == XPRS_LP_OPTIMAL || + xpress_status_ == XPRS_LP_UNFINISHED; + } } SolutionStatusProto XpressSolver::getLpSolutionStatus() const { - // TODO : put all statuses here switch (xpress_status_) { case XPRS_LP_OPTIMAL: + case XPRS_LP_UNFINISHED: return SOLUTION_STATUS_FEASIBLE; case XPRS_LP_INFEAS: + case XPRS_LP_CUTOFF: + case XPRS_LP_CUTOFF_IN_DUAL: + case XPRS_LP_NONCONVEX: return SOLUTION_STATUS_INFEASIBLE; + case XPRS_LP_UNSTARTED: case XPRS_LP_UNBOUNDED: + case XPRS_LP_UNSOLVED: return SOLUTION_STATUS_UNDETERMINED; default: return SOLUTION_STATUS_UNSPECIFIED; @@ -486,6 +556,19 @@ absl::StatusOr XpressSolver::GetSolveStats( SolveStatsProto solve_stats; CHECK_OK(util_time::EncodeGoogleApiProto(absl::Now() - start, solve_stats.mutable_solve_time())); + + // TODO : only run one of the following depending on algorithm? + // LP simplex iterations + { + ASSIGN_OR_RETURN(const int iters, xpress_->GetIntAttr(XPRS_SIMPLEXITER)); + solve_stats.set_simplex_iterations(iters); + } + // LP barrier iterations + { + ASSIGN_OR_RETURN(const int iters, xpress_->GetIntAttr(XPRS_BARITER)); + solve_stats.set_barrier_iterations(iters); + } + // TODO : complete these stats return solve_stats; } @@ -525,8 +608,10 @@ absl::StatusOr XpressSolver::ConvertTerminationReason( return CutoffTerminationProto( is_maximize_, "Objective worse than cutoff (XPRS_LP_CUTOFF)"); case XPRS_LP_UNFINISHED: - return LimitTerminationProto( - is_maximize_, LIMIT_UNSPECIFIED, best_primal_bound, best_dual_bound, + // TODO: add support for more limit types here (this only works for LP + // iterations limit for now) + return FeasibleTerminationProto( + is_maximize_, LIMIT_ITERATION, best_primal_bound, best_dual_bound, "Solve did not finish (XPRS_LP_UNFINISHED)"); case XPRS_LP_UNBOUNDED: return UnboundedTerminationProto(is_maximize_, diff --git a/ortools/math_opt/solvers/xpress_solver.h b/ortools/math_opt/solvers/xpress_solver.h index 62899ecad3e..4357effa5c4 100644 --- a/ortools/math_opt/solvers/xpress_solver.h +++ b/ortools/math_opt/solvers/xpress_solver.h @@ -57,7 +57,7 @@ class XpressSolver : public SolverInterface { MessageCallback message_cb, const CallbackRegistrationProto& callback_registration, Callback cb, const SolveInterrupter* interrupter) override; - absl::Status CallXpressSolve(bool enableOutput); + absl::Status CallXpressSolve(const SolveParametersProto& parameters); // Updates the problem (not implemented yet) absl::StatusOr Update(const ModelUpdateProto& model_update) override; diff --git a/ortools/math_opt/solvers/xpress_solver_test.cc b/ortools/math_opt/solvers/xpress_solver_test.cc index f038e3215cf..97094d3a21b 100644 --- a/ortools/math_opt/solvers/xpress_solver_test.cc +++ b/ortools/math_opt/solvers/xpress_solver_test.cc @@ -48,6 +48,28 @@ INSTANTIATE_TEST_SUITE_P( /*ensures_primal_ray=*/false, /*ensures_dual_ray=*/false, /*disallows_infeasible_or_unbounded=*/true))); +// TODO: implement missing LP features +INSTANTIATE_TEST_SUITE_P(XpressLpModelSolveParametersTest, + LpModelSolveParametersTest, + testing::Values(LpModelSolveParametersTestParameters( + SolverType::kXpress, /*exact_zeros=*/true, + /*supports_duals=*/true, + /*supports_primal_only_warm_starts=*/false))); + +// TODO: implement missing LP features +INSTANTIATE_TEST_SUITE_P( + GlopLpParameterTest, LpParameterTest, + testing::Values(LpParameterTestParams(SolverType::kXpress, + /*supports_simplex=*/true, + /*supports_barrier=*/true, + /*supports_first_order=*/false, + /*supports_random_seed=*/false, + /*supports_presolve=*/false, + /*supports_cutoff=*/false, + /*supports_objective_limit=*/false, + /*supports_best_bound_limit=*/false, + /*reports_limits=*/false))); + // TODO: implement message callbacks INSTANTIATE_TEST_SUITE_P(XpressMessageCallbackTest, MessageCallbackTest, testing::Values(MessageCallbackTestParams( @@ -97,7 +119,17 @@ INSTANTIATE_TEST_SUITE_P(XpressInfeasibleSubsystemTest, InfeasibleSubsystemTest, testing::Values(InfeasibleSubsystemTestParameters( {.solver_type = SolverType::kXpress, .support_menu = { - .supports_infeasible_subsystems = false}}))); + .supports_infeasible_subsystems = false}}))); + + +// TODO: When XPRESS supports MIP hints, enable the following test +GTEST_ALLOW_UNINSTANTIATED_PARAMETERIZED_TEST(MipSolutionHintTest); + +// TODO: When XPRESS supports MIP branch priorities, enable the following test +GTEST_ALLOW_UNINSTANTIATED_PARAMETERIZED_TEST(BranchPrioritiesTest); + +// TODO: When XPRESS supports lazy constraints, enable the following test +GTEST_ALLOW_UNINSTANTIATED_PARAMETERIZED_TEST(LazyConstraintsTest); } // namespace } // namespace math_opt diff --git a/ortools/xpress/environment.h b/ortools/xpress/environment.h index ae6e34a7c1a..ab5cc400d6a 100644 --- a/ortools/xpress/environment.h +++ b/ortools/xpress/environment.h @@ -420,6 +420,7 @@ absl::Status LoadXpressDynamicLibrary(std::string& xpresspath); #define XPRS_OBJSENSE 2008 #define XPRS_ROWS 1001 #define XPRS_SIMPLEXITER 1009 +#define XPRS_BARITER 5001 #define XPRS_LPSTATUS 1010 #define XPRS_MIPSTATUS 1011 #define XPRS_NODES 1013 From 38a13bacf2c8083be00eda6f862acb4bcfe8225b Mon Sep 17 00:00:00 2001 From: Peter Mitri Date: Fri, 6 Dec 2024 17:26:43 +0100 Subject: [PATCH 11/14] Enable all tests required by Google. Improve some features and clean up the code. --- ortools/math_opt/solvers/xpress/g_xpress.cc | 26 ++- ortools/math_opt/solvers/xpress/g_xpress.h | 5 +- ortools/math_opt/solvers/xpress_solver.cc | 195 +++++++++--------- ortools/math_opt/solvers/xpress_solver.h | 4 +- .../math_opt/solvers/xpress_solver_test.cc | 142 +++++++++++-- 5 files changed, 250 insertions(+), 122 deletions(-) diff --git a/ortools/math_opt/solvers/xpress/g_xpress.cc b/ortools/math_opt/solvers/xpress/g_xpress.cc index 00d7baf0390..d7bf76ec44c 100644 --- a/ortools/math_opt/solvers/xpress/g_xpress.cc +++ b/ortools/math_opt/solvers/xpress/g_xpress.cc @@ -32,8 +32,6 @@ #include "ortools/base/status_macros.h" #include "ortools/xpress/environment.h" -// TODO : try to remove all CHECK_EQ and replace them with absl Status - namespace operations_research::math_opt { constexpr int kXpressOk = 0; @@ -88,18 +86,17 @@ absl::Status Xpress::AddVars(const absl::Span vbegin, const absl::Span lb, const absl::Span ub, const absl::Span vtype) { - CHECK_EQ(vind.size(), vval.size()); const int num_vars = static_cast(lb.size()); - CHECK_EQ(ub.size(), num_vars); - CHECK_EQ(vtype.size(), num_vars); + if (vind.size() != vval.size() || ub.size() != num_vars || + vtype.size() != num_vars || (!obj.empty() && obj.size() != num_vars) || + (!vbegin.empty() && vbegin.size() != num_vars)) { + return absl::InvalidArgumentError( + "Xpress::AddVars arguments are of inconsistent sizes"); + } double* c_obj = nullptr; if (!obj.empty()) { - CHECK_EQ(obj.size(), num_vars); c_obj = const_cast(obj.data()); } - if (!vbegin.empty()) { - CHECK_EQ(vbegin.size(), num_vars); - } // TODO: look into int64 support for number of vars (use XPRSaddcols64) return ToStatus(XPRSaddcols(xpress_model_, num_vars, 0, c_obj, nullptr, nullptr, nullptr, lb.data(), ub.data())); @@ -109,7 +106,10 @@ absl::Status Xpress::AddConstrs(const absl::Span sense, const absl::Span rhs, const absl::Span rng) { const int num_cons = static_cast(sense.size()); - CHECK_EQ(rhs.size(), num_cons); + if (rhs.size() != num_cons) { + return absl::InvalidArgumentError( + "RHS must have one element per constraint."); + } return ToStatus(XPRSaddrows(xpress_model_, num_cons, 0, sense.data(), rhs.data(), rng.data(), NULL, NULL, NULL)); } @@ -235,6 +235,12 @@ absl::Status Xpress::GetBasis(std::vector& rowBasis, XPRSgetbasis(xpress_model_, rowBasis.data(), colBasis.data())); } +absl::Status Xpress::SetStartingBasis(std::vector& rowBasis, + std::vector& colBasis) const { + return ToStatus( + XPRSloadbasis(xpress_model_, rowBasis.data(), colBasis.data())); +} + absl::StatusOr> Xpress::GetVarLb() const { int nVars = GetNumberOfVariables(); std::vector bounds; diff --git a/ortools/math_opt/solvers/xpress/g_xpress.h b/ortools/math_opt/solvers/xpress/g_xpress.h index 83b5cddb6c6..a18ebabd1d2 100644 --- a/ortools/math_opt/solvers/xpress/g_xpress.h +++ b/ortools/math_opt/solvers/xpress/g_xpress.h @@ -83,7 +83,10 @@ class Xpress { absl::StatusOr> GetPrimalValues() const; absl::StatusOr> GetConstraintDuals() const; absl::StatusOr> GetReducedCostValues() const; - absl::Status GetBasis(std::vector& rowBasis, std::vector& colBasis) const; + absl::Status GetBasis(std::vector& rowBasis, + std::vector& colBasis) const; + absl::Status SetStartingBasis(std::vector& rowBasis, + std::vector& colBasis) const; static void XPRS_CC printXpressMessage(XPRSprob prob, void* data, const char* sMsg, int nLen, diff --git a/ortools/math_opt/solvers/xpress_solver.cc b/ortools/math_opt/solvers/xpress_solver.cc index 6f11b0633db..7538841b542 100644 --- a/ortools/math_opt/solvers/xpress_solver.cc +++ b/ortools/math_opt/solvers/xpress_solver.cc @@ -18,6 +18,7 @@ #include "ortools/base/protoutil.h" #include "ortools/base/status_macros.h" #include "ortools/math_opt/core/math_opt_proto_utils.h" +#include "ortools/math_opt/core/sparse_vector_view.h" #include "ortools/math_opt/cpp/solve_result.h" #include "ortools/math_opt/validators/callback_validator.h" #include "ortools/port/proto_utils.h" @@ -26,7 +27,7 @@ namespace operations_research { namespace math_opt { namespace { -// TODO: add all unsupported parameter combinations here + absl::Status CheckParameters(const SolveParametersProto& parameters) { std::vector warnings; if (parameters.has_threads() && parameters.threads() > 1) { @@ -59,14 +60,14 @@ absl::Status CheckParameters(const SolveParametersProto& parameters) { } // namespace constexpr SupportedProblemStructures kXpressSupportedStructures = { - .integer_variables = SupportType::kNotImplemented, - .multi_objectives = SupportType::kNotImplemented, - .quadratic_objectives = SupportType::kNotImplemented, - .quadratic_constraints = SupportType::kNotImplemented, - .second_order_cone_constraints = SupportType::kNotImplemented, - .sos1_constraints = SupportType::kNotImplemented, - .sos2_constraints = SupportType::kNotImplemented, - .indicator_constraints = SupportType::kNotImplemented}; + .integer_variables = SupportType::kNotSupported, + .multi_objectives = SupportType::kNotSupported, + .quadratic_objectives = SupportType::kNotSupported, + .quadratic_constraints = SupportType::kNotSupported, + .second_order_cone_constraints = SupportType::kNotSupported, + .sos1_constraints = SupportType::kNotSupported, + .sos2_constraints = SupportType::kNotSupported, + .indicator_constraints = SupportType::kNotSupported}; absl::StatusOr> XpressSolver::New( const ModelProto& input_model, const SolverInterface::InitArgs& init_args) { @@ -118,7 +119,7 @@ absl::Status XpressSolver::AddNewVariables( variable_type)); // Not adding names for performance (have to call XPRSaddnames) - // TODO : keep names in a cache and add them when needed? + // TODO: keep names in a cache and add them when needed? return absl::OkStatus(); } @@ -175,13 +176,6 @@ absl::Status XpressSolver::AddNewLinearConstraints( } absl::Status XpressSolver::AddSingleObjective(const ObjectiveProto& objective) { - // TODO: reactivate the following code after figuring out why the condition is - // true for LP tests - /*if (objective.has_quadratic_coefficients()) { - return absl::UnimplementedError( - "Quadratic objectives are not yet implemented in XPRESS solver " - "interface."); - }*/ std::vector index; index.reserve(objective.linear_coefficients().ids_size()); for (const int64_t id : objective.linear_coefficients().ids()) { @@ -215,15 +209,10 @@ absl::StatusOr XpressSolver::Solve( MessageCallback message_cb, const CallbackRegistrationProto& callback_registration, Callback cb, const SolveInterrupter* interrupter) { + RETURN_IF_ERROR(ModelSolveParametersAreSupported( + model_parameters, kXpressSupportedStructures, "XPRESS")); const absl::Time start = absl::Now(); - // TODO: set solve parameters - // TODO: set basis - // TODO: set hints - // TODO: set branching properties - // TODO: set lazy constraints - // TODO: add interrupter using xpress_->Terminate(); - RETURN_IF_ERROR(CheckRegisteredCallbackEvents(callback_registration, /*supported_events=*/{})); @@ -236,6 +225,11 @@ absl::StatusOr XpressSolver::Solve( RETURN_IF_ERROR(inv_bounds.ToStatus()); } + // Set initial basis + if (model_parameters.has_initial_basis()) { + RETURN_IF_ERROR(SetXpressStartingBasis(model_parameters.initial_basis())); + } + RETURN_IF_ERROR(CallXpressSolve(parameters)) << "Error during XPRESS solve"; ASSIGN_OR_RETURN(SolveResultProto solve_result, @@ -244,6 +238,20 @@ absl::StatusOr XpressSolver::Solve( return solve_result; } +inline std::string GetOptimizationFlags( + const SolveParametersProto& parameters) { + switch (parameters.lp_algorithm()) { + case LP_ALGORITHM_PRIMAL_SIMPLEX: + return "p"; + case LP_ALGORITHM_DUAL_SIMPLEX: + return "d"; + case LP_ALGORITHM_BARRIER: + return "b"; + default: + // this makes XPRESS use default algorithm (XPRS_DEFAULTALG) + return ""; + } +} absl::Status XpressSolver::CallXpressSolve( const SolveParametersProto& parameters) { // Enable screen output right before solve @@ -255,50 +263,10 @@ absl::Status XpressSolver::CallXpressSolve( if (is_mip_) { ASSIGN_OR_RETURN(xpress_status_, xpress_->MipOptimizeAndGetStatus()); } else { - std::string flags; - // TODO: put "p", "d", "b" in environment.h? - // TODO: refactor this, it is ugly - switch (parameters.lp_algorithm()) { - case LP_ALGORITHM_PRIMAL_SIMPLEX: - flags = "p"; - if (parameters.has_iteration_limit()) { - RETURN_IF_ERROR(xpress_->SetIntAttr(XPRS_LPITERLIMIT, - parameters.iteration_limit())); - } - break; - case LP_ALGORITHM_DUAL_SIMPLEX: - flags = "d"; - if (parameters.has_iteration_limit()) { - RETURN_IF_ERROR(xpress_->SetIntAttr(XPRS_LPITERLIMIT, - parameters.iteration_limit())); - } - break; - case LP_ALGORITHM_BARRIER: - flags = "b"; - if (parameters.has_iteration_limit()) { - RETURN_IF_ERROR(xpress_->SetIntAttr(XPRS_BARITERLIMIT, - parameters.iteration_limit())); - // we also have to set simplex iter limits because barrier iterations - // call simplex iterations - // TODO: verify this - RETURN_IF_ERROR(xpress_->SetIntAttr(XPRS_LPITERLIMIT, - parameters.iteration_limit())); - } - break; - default: - // this makes XPRESS use default algorithm (XPRS_DEFAULTALG) - flags = ""; - if (parameters.has_iteration_limit()) { - // TODO: fetch the default alg and set its limit? or default to primal - // simplex - RETURN_IF_ERROR(xpress_->SetIntAttr(XPRS_LPITERLIMIT, - parameters.iteration_limit())); - RETURN_IF_ERROR(xpress_->SetIntAttr(XPRS_BARITERLIMIT, - parameters.iteration_limit())); - } - break; - } - ASSIGN_OR_RETURN(xpress_status_, xpress_->LpOptimizeAndGetStatus(flags)); + RETURN_IF_ERROR(SetLpIterLimits(parameters)) + << "Could not set iteration limits."; + ASSIGN_OR_RETURN(xpress_status_, xpress_->LpOptimizeAndGetStatus( + GetOptimizationFlags(parameters))); } // Post-solve if (!(is_mip_ ? (xpress_status_ == XPRS_MIP_OPTIMAL) @@ -313,6 +281,24 @@ absl::Status XpressSolver::CallXpressSolve( return absl::OkStatus(); } +absl::Status XpressSolver::SetLpIterLimits( + const SolveParametersProto& parameters) { + // If the user has set no limits, we still have to reset the limits + // explicitly to their default values, else the parameters could be kept + // between solves. XPRESS default value for XPRS_LPITERLIMIT is 2147483647 + int lp_iter_limit = parameters.has_iteration_limit() + ? parameters.iteration_limit() + : 2147483647; + // XPRESS default value for XPRS_BARITERLIMIT is 500 + int bar_iter_limit = + parameters.has_iteration_limit() ? parameters.iteration_limit() : 500; + RETURN_IF_ERROR(xpress_->SetIntAttr(XPRS_LPITERLIMIT, lp_iter_limit)) + << "Could not set XPRS_LPITERLIMIT"; + RETURN_IF_ERROR(xpress_->SetIntAttr(XPRS_BARITERLIMIT, bar_iter_limit)) + << "Could not set XPRS_BARITERLIMIT"; + return absl::OkStatus(); +} + absl::StatusOr XpressSolver::ExtractSolveResultProto( absl::Time start, const ModelSolveParametersProto& model_parameters) { SolveResultProto result; @@ -396,6 +382,9 @@ SolutionStatusProto XpressSolver::getLpSolutionStatus() const { switch (xpress_status_) { case XPRS_LP_OPTIMAL: case XPRS_LP_UNFINISHED: + // TODO: XPRESS returns XPRS_LP_UNFINISHED even if it found no solution + // maybe we should adapt this code (try it out with LpIncompleteSolveTest, + // by setting "initial basis" support flags to true) return SOLUTION_STATUS_FEASIBLE; case XPRS_LP_INFEAS: case XPRS_LP_CUTOFF: @@ -473,14 +462,19 @@ XpressSolver::GetConvexDualSolutionIfAvailable( .feasible_solution_exists = dual_feasible_solution_exists}; } -inline BasisStatusProto ConvertVariableStatus(const int status) { +inline BasisStatusProto XpressToMathOptBasisStatus(const int status, + bool isConstraint) { + // XPRESS row basis status is that of the slack variable + // For example, if the slack variable is at LB, the constraint is at UB switch (status) { case XPRS_BASIC: return BASIS_STATUS_BASIC; case XPRS_AT_LOWER: - return BASIS_STATUS_AT_LOWER_BOUND; + return isConstraint ? BASIS_STATUS_AT_UPPER_BOUND + : BASIS_STATUS_AT_LOWER_BOUND; case XPRS_AT_UPPER: - return BASIS_STATUS_AT_UPPER_BOUND; + return isConstraint ? BASIS_STATUS_AT_LOWER_BOUND + : BASIS_STATUS_AT_UPPER_BOUND; case XPRS_FREE_SUPER: return BASIS_STATUS_FREE; default: @@ -488,23 +482,41 @@ inline BasisStatusProto ConvertVariableStatus(const int status) { } } -inline BasisStatusProto ConvertConstraintStatus(const int status) { +inline int MathOptToXpressBasisStatus(const BasisStatusProto status, + bool isConstraint) { // XPRESS row basis status is that of the slack variable // For example, if the slack variable is at LB, the constraint is at UB switch (status) { - case XPRS_BASIC: - return BASIS_STATUS_BASIC; - case XPRS_AT_LOWER: - return BASIS_STATUS_AT_UPPER_BOUND; - case XPRS_AT_UPPER: - return BASIS_STATUS_AT_LOWER_BOUND; - case XPRS_FREE_SUPER: - return BASIS_STATUS_FREE; + case BASIS_STATUS_BASIC: + return XPRS_BASIC; + case BASIS_STATUS_AT_LOWER_BOUND: + return isConstraint ? XPRS_AT_UPPER : XPRS_AT_LOWER; + case BASIS_STATUS_AT_UPPER_BOUND: + return isConstraint ? XPRS_AT_LOWER : XPRS_AT_UPPER; + case BASIS_STATUS_FREE: + return XPRS_FREE_SUPER; default: - return BASIS_STATUS_UNSPECIFIED; + return XPRS_FREE_SUPER; } } +absl::Status XpressSolver::SetXpressStartingBasis(const BasisProto& basis) { + std::vector xpress_var_basis_status(xpress_->GetNumberOfVariables()); + for (const auto [id, value] : MakeView(basis.variable_status())) { + xpress_var_basis_status[variables_map_.at(id)] = + MathOptToXpressBasisStatus(static_cast(value), false); + } + std::vector xpress_constr_basis_status( + xpress_->GetNumberOfConstraints()); + for (const auto [id, value] : MakeView(basis.constraint_status())) { + xpress_constr_basis_status[linear_constraints_map_.at(id) + .constraint_index] = + MathOptToXpressBasisStatus(static_cast(value), true); + } + return xpress_->SetStartingBasis(xpress_constr_basis_status, + xpress_var_basis_status); +} + absl::StatusOr> XpressSolver::GetBasisIfAvailable() { BasisProto basis; @@ -517,8 +529,8 @@ absl::StatusOr> XpressSolver::GetBasisIfAvailable() { // Variable basis for (auto [variable_id, xprs_variable_index] : variables_map_) { basis.mutable_variable_status()->add_ids(variable_id); - const BasisStatusProto variable_status = - ConvertVariableStatus(xprs_variable_basis_status[xprs_variable_index]); + const BasisStatusProto variable_status = XpressToMathOptBasisStatus( + xprs_variable_basis_status[xprs_variable_index], false); if (variable_status == BASIS_STATUS_UNSPECIFIED) { return absl::InternalError( absl::StrCat("Invalid Xpress variable basis status: ", @@ -530,8 +542,8 @@ absl::StatusOr> XpressSolver::GetBasisIfAvailable() { // Constraint basis for (auto [constraint_id, xprs_ct_index] : linear_constraints_map_) { basis.mutable_constraint_status()->add_ids(constraint_id); - const BasisStatusProto status = ConvertConstraintStatus( - xprs_constraint_basis_status[xprs_ct_index.constraint_index]); + const BasisStatusProto status = XpressToMathOptBasisStatus( + xprs_constraint_basis_status[xprs_ct_index.constraint_index], true); if (status == BASIS_STATUS_UNSPECIFIED) { return absl::InternalError(absl::StrCat( "Invalid Xpress constraint basis status: ", @@ -557,7 +569,6 @@ absl::StatusOr XpressSolver::GetSolveStats( CHECK_OK(util_time::EncodeGoogleApiProto(absl::Now() - start, solve_stats.mutable_solve_time())); - // TODO : only run one of the following depending on algorithm? // LP simplex iterations { ASSIGN_OR_RETURN(const int iters, xpress_->GetIntAttr(XPRS_SIMPLEXITER)); @@ -569,7 +580,7 @@ absl::StatusOr XpressSolver::GetSolveStats( solve_stats.set_barrier_iterations(iters); } - // TODO : complete these stats + // TODO: complete these stats return solve_stats; } @@ -639,10 +650,7 @@ absl::StatusOr XpressSolver::ConvertTerminationReason( absl::StatusOr XpressSolver::Update( const ModelUpdateProto& model_update) { - if (!UpdateIsSupported(model_update, kXpressSupportedStructures)) { - return false; - } - // TODO: implement Update + // Not implemented yet return false; } @@ -650,7 +658,6 @@ absl::StatusOr XpressSolver::ComputeInfeasibleSubsystem(const SolveParametersProto& parameters, MessageCallback message_cb, const SolveInterrupter* interrupter) { - // TODO: implement this return absl::UnimplementedError( "XpressSolver cannot compute infeasible subsystem yet"); } @@ -666,7 +673,9 @@ absl::StatusOr XpressSolver::ListInvertedBounds() const { } } } - // TODO: we better fetch constraint lb & ub from xpress_ to ensure sync + // We could have used XPRSgetrhsrange to check if one is negative. However, + // XPRSaddrows ignores the sign of the RHS range and takes the absolute value + // in all cases. So we need to do the following checks on the internal model. for (const auto& [id, cstr_data] : linear_constraints_map_) { if (cstr_data.lower_bound > cstr_data.upper_bound) { inverted_bounds.linear_constraints.push_back(id); diff --git a/ortools/math_opt/solvers/xpress_solver.h b/ortools/math_opt/solvers/xpress_solver.h index 4357effa5c4..0ce9cf31c94 100644 --- a/ortools/math_opt/solvers/xpress_solver.h +++ b/ortools/math_opt/solvers/xpress_solver.h @@ -193,10 +193,10 @@ class XpressSolver : public SolverInterface { int get_model_index(const LinearConstraintData& index) const { return index.constraint_index; } - SolutionStatusProto getLpSolutionStatus() const; - absl::StatusOr ListInvertedBounds() const; + absl::Status SetXpressStartingBasis(const BasisProto& basis); + absl::Status SetLpIterLimits(const SolveParametersProto& parameters); // Fields to track the number of Xpress variables and constraints. These // quantities are updated immediately after adding or removing to the model, diff --git a/ortools/math_opt/solvers/xpress_solver_test.cc b/ortools/math_opt/solvers/xpress_solver_test.cc index 97094d3a21b..64798bed72c 100644 --- a/ortools/math_opt/solvers/xpress_solver_test.cc +++ b/ortools/math_opt/solvers/xpress_solver_test.cc @@ -39,7 +39,6 @@ namespace math_opt { namespace { using testing::ValuesIn; -// TODO: implement missing LP features INSTANTIATE_TEST_SUITE_P( XpressSolverLpTest, SimpleLpTest, testing::Values(SimpleLpTestParameters( @@ -48,7 +47,6 @@ INSTANTIATE_TEST_SUITE_P( /*ensures_primal_ray=*/false, /*ensures_dual_ray=*/false, /*disallows_infeasible_or_unbounded=*/true))); -// TODO: implement missing LP features INSTANTIATE_TEST_SUITE_P(XpressLpModelSolveParametersTest, LpModelSolveParametersTest, testing::Values(LpModelSolveParametersTestParameters( @@ -56,9 +54,8 @@ INSTANTIATE_TEST_SUITE_P(XpressLpModelSolveParametersTest, /*supports_duals=*/true, /*supports_primal_only_warm_starts=*/false))); -// TODO: implement missing LP features INSTANTIATE_TEST_SUITE_P( - GlopLpParameterTest, LpParameterTest, + XpressLpParameterTest, LpParameterTest, testing::Values(LpParameterTestParams(SolverType::kXpress, /*supports_simplex=*/true, /*supports_barrier=*/true, @@ -70,7 +67,29 @@ INSTANTIATE_TEST_SUITE_P( /*supports_best_bound_limit=*/false, /*reports_limits=*/false))); -// TODO: implement message callbacks +INSTANTIATE_TEST_SUITE_P( + XpressPrimalSimplexLpIncompleteSolveTest, LpIncompleteSolveTest, + testing::Values(LpIncompleteSolveTestParams( + SolverType::kXpress, + /*lp_algorithm=*/LPAlgorithm::kPrimalSimplex, + /*supports_iteration_limit=*/true, /*supports_initial_basis=*/false, + /*supports_incremental_solve=*/false, /*supports_basis=*/true, + /*supports_presolve=*/false, /*check_primal_objective=*/true, + /*primal_solution_status_always_set=*/true, + /*dual_solution_status_always_set=*/true))); +INSTANTIATE_TEST_SUITE_P( + XpressDualSimplexLpIncompleteSolveTest, LpIncompleteSolveTest, + testing::Values(LpIncompleteSolveTestParams( + SolverType::kXpress, + /*lp_algorithm=*/LPAlgorithm::kDualSimplex, + /*supports_iteration_limit=*/true, /*supports_initial_basis=*/false, + /*supports_incremental_solve=*/false, /*supports_basis=*/true, + /*supports_presolve=*/false, /*check_primal_objective=*/true, + /*primal_solution_status_always_set=*/true, + /*dual_solution_status_always_set=*/true))); + +GTEST_ALLOW_UNINSTANTIATED_PARAMETERIZED_TEST(IncrementalLpTest); + INSTANTIATE_TEST_SUITE_P(XpressMessageCallbackTest, MessageCallbackTest, testing::Values(MessageCallbackTestParams( SolverType::kXpress, @@ -78,7 +97,6 @@ INSTANTIATE_TEST_SUITE_P(XpressMessageCallbackTest, MessageCallbackTest, /*support_interrupter=*/false, /*integer_variables=*/false, ""))); -// TODO: implement callbacks INSTANTIATE_TEST_SUITE_P( XpressCallbackTest, CallbackTest, testing::Values(CallbackTestParams(SolverType::kXpress, @@ -101,7 +119,6 @@ InvalidParameterTestParams InvalidThreadsParameters() { {"only supports parameters.threads = 1"}); } -// TODO: add all invalid parameters combinations INSTANTIATE_TEST_SUITE_P(XpressInvalidParameterTest, InvalidParameterTest, ValuesIn({InvalidThreadsParameters()})); @@ -111,26 +128,119 @@ INSTANTIATE_TEST_SUITE_P(XpressGenericTest, GenericTest, /*integer_variables=*/false, /*expected_log=*/"Optimal solution found"))); -// TODO: When XPRESS callbacks are supported, enable this test. GTEST_ALLOW_UNINSTANTIATED_PARAMETERIZED_TEST(TimeLimitTest); -// TODO: implement IIS support INSTANTIATE_TEST_SUITE_P(XpressInfeasibleSubsystemTest, InfeasibleSubsystemTest, testing::Values(InfeasibleSubsystemTestParameters( {.solver_type = SolverType::kXpress, .support_menu = { - .supports_infeasible_subsystems = false}}))); + .supports_infeasible_subsystems = false}}))); - -// TODO: When XPRESS supports MIP hints, enable the following test +GTEST_ALLOW_UNINSTANTIATED_PARAMETERIZED_TEST(IpModelSolveParametersTest); +GTEST_ALLOW_UNINSTANTIATED_PARAMETERIZED_TEST(IpParameterTest); +GTEST_ALLOW_UNINSTANTIATED_PARAMETERIZED_TEST(LargeInstanceIpParameterTest); +GTEST_ALLOW_UNINSTANTIATED_PARAMETERIZED_TEST(SimpleMipTest); +GTEST_ALLOW_UNINSTANTIATED_PARAMETERIZED_TEST(IncrementalMipTest); GTEST_ALLOW_UNINSTANTIATED_PARAMETERIZED_TEST(MipSolutionHintTest); - -// TODO: When XPRESS supports MIP branch priorities, enable the following test GTEST_ALLOW_UNINSTANTIATED_PARAMETERIZED_TEST(BranchPrioritiesTest); - -// TODO: When XPRESS supports lazy constraints, enable the following test GTEST_ALLOW_UNINSTANTIATED_PARAMETERIZED_TEST(LazyConstraintsTest); +LogicalConstraintTestParameters GetXpressLogicalConstraintTestParameters() { + return LogicalConstraintTestParameters( + SolverType::kXpress, SolveParameters(), + /*supports_integer_variables=*/false, + /*supports_sos1=*/false, + /*supports_sos2=*/false, + /*supports_indicator_constraints=*/false, + /*supports_incremental_add_and_deletes=*/false, + /*supports_incremental_variable_deletions=*/false, + /*supports_deleting_indicator_variables=*/false, + /*supports_updating_binary_variables=*/false); +} + +INSTANTIATE_TEST_SUITE_P( + XpressSimpleLogicalConstraintTest, SimpleLogicalConstraintTest, + testing::Values(GetXpressLogicalConstraintTestParameters())); +INSTANTIATE_TEST_SUITE_P( + XpressIncrementalLogicalConstraintTest, IncrementalLogicalConstraintTest, + testing::Values(GetXpressLogicalConstraintTestParameters())); + +MultiObjectiveTestParameters GetXpressMultiObjectiveTestParameters() { + return MultiObjectiveTestParameters( + /*solver_type=*/SolverType::kXpress, /*parameters=*/SolveParameters(), + /*supports_auxiliary_objectives=*/false, + /*supports_incremental_objective_add_and_delete=*/false, + /*supports_incremental_objective_modification=*/false, + /*supports_integer_variables=*/false); +} + +INSTANTIATE_TEST_SUITE_P( + XpressSimpleMultiObjectiveTest, SimpleMultiObjectiveTest, + testing::Values(GetXpressMultiObjectiveTestParameters())); + +INSTANTIATE_TEST_SUITE_P( + XpressIncrementalMultiObjectiveTest, IncrementalMultiObjectiveTest, + testing::Values(GetXpressMultiObjectiveTestParameters())); + +QpTestParameters GetXpressQpTestParameters() { + return QpTestParameters(SolverType::kXpress, SolveParameters(), + /*qp_support=*/QpSupportType::kNoQpSupport, + /*supports_incrementalism_not_modifying_qp=*/false, + /*supports_qp_incrementalism=*/false, + /*use_integer_variables=*/false); +} +INSTANTIATE_TEST_SUITE_P(XpressSimpleQpTest, SimpleQpTest, + testing::Values(GetXpressQpTestParameters())); +INSTANTIATE_TEST_SUITE_P(XpressIncrementalQpTest, IncrementalQpTest, + testing::Values(GetXpressQpTestParameters())); +GTEST_ALLOW_UNINSTANTIATED_PARAMETERIZED_TEST(QpDualsTest); + +QcTestParameters GetXpressQcTestParameters() { + return QcTestParameters(SolverType::kXpress, SolveParameters(), + /*supports_qc=*/false, + /*supports_incremental_add_and_deletes=*/false, + /*supports_incremental_variable_deletions=*/false, + /*use_integer_variables=*/false); +} +INSTANTIATE_TEST_SUITE_P(XpressSimpleQcTest, SimpleQcTest, + testing::Values(GetXpressQcTestParameters())); +INSTANTIATE_TEST_SUITE_P(XpressIncrementalQcTest, IncrementalQcTest, + testing::Values(GetXpressQcTestParameters())); +GTEST_ALLOW_UNINSTANTIATED_PARAMETERIZED_TEST(QcDualsTest); + +SecondOrderConeTestParameters GetXpressSecondOrderConeTestParameters() { + return SecondOrderConeTestParameters( + SolverType::kXpress, SolveParameters(), + /*supports_soc_constraints=*/false, + /*supports_incremental_add_and_deletes=*/false); +} +INSTANTIATE_TEST_SUITE_P( + XpressSimpleSecondOrderConeTest, SimpleSecondOrderConeTest, + testing::Values(GetXpressSecondOrderConeTestParameters())); +INSTANTIATE_TEST_SUITE_P( + XpressIncrementalSecondOrderConeTest, IncrementalSecondOrderConeTest, + testing::Values(GetXpressSecondOrderConeTestParameters())); + +std::vector MakeStatusTestConfigs() { + std::vector test_parameters; + for (const auto algorithm : std::vector>( + {std::nullopt, LPAlgorithm::kBarrier, LPAlgorithm::kPrimalSimplex, + LPAlgorithm::kDualSimplex})) { + SolveParameters solve_parameters = {.lp_algorithm = algorithm}; + test_parameters.push_back(StatusTestParameters( + SolverType::kXpress, solve_parameters, + /*disallow_primal_or_dual_infeasible=*/false, + /*supports_iteration_limit=*/false, + /*use_integer_variables=*/false, + /*supports_node_limit=*/false, + /*support_interrupter=*/false, /*supports_one_thread=*/true)); + } + return test_parameters; +} + +INSTANTIATE_TEST_SUITE_P(XpressStatusTest, StatusTest, + ValuesIn(MakeStatusTestConfigs())); + } // namespace } // namespace math_opt } // namespace operations_research \ No newline at end of file From 1fb56d822e4b813dbb2246f0eb6659948d9bb10e Mon Sep 17 00:00:00 2001 From: Peter Mitri Date: Fri, 3 Jan 2025 16:50:18 +0100 Subject: [PATCH 12/14] review from Daniel --- ortools/math_opt/solvers/xpress/g_xpress.cc | 96 +++++++++++++++------ ortools/math_opt/solvers/xpress/g_xpress.h | 9 +- ortools/math_opt/solvers/xpress_solver.cc | 32 ++++--- ortools/xpress/environment.cc | 10 +++ ortools/xpress/environment.h | 6 ++ 5 files changed, 110 insertions(+), 43 deletions(-) diff --git a/ortools/math_opt/solvers/xpress/g_xpress.cc b/ortools/math_opt/solvers/xpress/g_xpress.cc index d7bf76ec44c..1488de895de 100644 --- a/ortools/math_opt/solvers/xpress/g_xpress.cc +++ b/ortools/math_opt/solvers/xpress/g_xpress.cc @@ -33,6 +33,11 @@ #include "ortools/xpress/environment.h" namespace operations_research::math_opt { + +namespace { +bool checkInt32Overflow(const unsigned int value) { return value > INT32_MAX; } +} // namespace + constexpr int kXpressOk = 0; absl::Status Xpress::ToStatus(const int xprs_err, @@ -41,22 +46,28 @@ absl::Status Xpress::ToStatus(const int xprs_err, return absl::OkStatus(); } char errmsg[512]; - XPRSgetlasterror(xpress_model_, errmsg); - return util::StatusBuilder(code) - << "Xpress error code: " << xprs_err << ", message: " << errmsg; + int status = XPRSgetlasterror(xpress_model_, errmsg); + if (status == kXpressOk) { + return util::StatusBuilder(code) + << "Xpress error code: " << xprs_err << ", message: " << errmsg; + } + return util::StatusBuilder(code) << "Xpress error code: " << xprs_err + << " (message could not be fetched)"; } -Xpress::Xpress(XPRSprob& model) : xpress_model_(ABSL_DIE_IF_NULL(model)) {} +Xpress::Xpress(XPRSprob& model) : xpress_model_(ABSL_DIE_IF_NULL(model)) { + initIntControlDefaults(); +} absl::StatusOr> Xpress::New( const std::string& model_name) { bool correctlyLoaded = initXpressEnv(); CHECK(correctlyLoaded); - XPRSprob* model; - CHECK_EQ(kXpressOk, XPRScreateprob(model)); + XPRSprob model; + CHECK_EQ(kXpressOk, XPRScreateprob(&model)); DCHECK(model != nullptr); // should not be NULL if status=0 - CHECK_EQ(kXpressOk, XPRSaddcbmessage(*model, printXpressMessage, NULL, 0)); - return absl::WrapUnique(new Xpress(*model)); + CHECK_EQ(kXpressOk, XPRSaddcbmessage(model, printXpressMessage, nullptr, 0)); + return absl::WrapUnique(new Xpress(model)); } void XPRS_CC Xpress::printXpressMessage(XPRSprob prob, void* data, @@ -72,6 +83,13 @@ Xpress::~Xpress() { CHECK_EQ(kXpressOk, XPRSfree()); } +void Xpress::initIntControlDefaults() { + std::vector controls = {XPRS_LPITERLIMIT, XPRS_BARITERLIMIT}; + for (auto control : controls) { + int_control_defaults_[control] = GetIntControl(control).value(); + } +} + absl::Status Xpress::AddVars(const absl::Span obj, const absl::Span lb, const absl::Span ub, @@ -86,13 +104,17 @@ absl::Status Xpress::AddVars(const absl::Span vbegin, const absl::Span lb, const absl::Span ub, const absl::Span vtype) { + if (checkInt32Overflow(lb.size())) { + return absl::InvalidArgumentError( + "XPRESS cannot handle more than 2^31 variables"); + } const int num_vars = static_cast(lb.size()); if (vind.size() != vval.size() || ub.size() != num_vars || vtype.size() != num_vars || (!obj.empty() && obj.size() != num_vars) || (!vbegin.empty() && vbegin.size() != num_vars)) { - return absl::InvalidArgumentError( - "Xpress::AddVars arguments are of inconsistent sizes"); - } + return absl::InvalidArgumentError( + "Xpress::AddVars arguments are of inconsistent sizes"); + } double* c_obj = nullptr; if (!obj.empty()) { c_obj = const_cast(obj.data()); @@ -127,20 +149,16 @@ absl::Status Xpress::SetObjective(bool maximize, double offset, << "Failed to set objective offset in XPRESS"; const int n_cols = static_cast(colind.size()); - auto c_colind = const_cast(colind.data()); - auto c_values = const_cast(values.data()); - return ToStatus(XPRSchgobj(xpress_model_, n_cols, c_colind, c_values)); + return ToStatus( + XPRSchgobj(xpress_model_, n_cols, colind.data(), values.data())); } absl::Status Xpress::ChgCoeffs(absl::Span rowind, absl::Span colind, absl::Span values) { - const int n_coefs = static_cast(rowind.size()); - auto c_rowind = const_cast(rowind.data()); - auto c_colind = const_cast(colind.data()); - auto c_values = const_cast(values.data()); - return ToStatus( - XPRSchgmcoef(xpress_model_, n_coefs, c_rowind, c_colind, c_values)); + const long n_coefs = static_cast(rowind.size()); + return ToStatus(XPRSchgmcoef64(xpress_model_, n_coefs, rowind.data(), + colind.data(), values.data())); } absl::StatusOr Xpress::LpOptimizeAndGetStatus(std::string flags) { @@ -168,6 +186,28 @@ absl::StatusOr Xpress::MipOptimizeAndGetStatus() { void Xpress::Terminate() { XPRSinterrupt(xpress_model_, XPRS_STOP_USER); }; +absl::StatusOr Xpress::GetIntControl(int control) const { + int result; + RETURN_IF_ERROR(ToStatus(XPRSgetintcontrol(xpress_model_, control, &result))) + << "Error getting Xpress int control: " << control; + return result; +} + +absl::Status Xpress::SetIntControl(int control, int value) { + return ToStatus(XPRSsetintcontrol(xpress_model_, control, value)); +} + +absl::Status Xpress::ResetIntControl(int control) { + if (int_control_defaults_.count(control)) { + return ToStatus(XPRSsetintcontrol(xpress_model_, control, + int_control_defaults_[control])); + } + return absl::InvalidArgumentError( + "Default value unknown for control " + std::to_string(control) + + ", consider adding it to Xpress::initIntControlDefaults"); +} + + absl::StatusOr Xpress::GetIntAttr(int attribute) const { int result; RETURN_IF_ERROR(ToStatus(XPRSgetintattrib(xpress_model_, attribute, &result))) @@ -175,10 +215,6 @@ absl::StatusOr Xpress::GetIntAttr(int attribute) const { return result; } -absl::Status Xpress::SetIntAttr(int attribute, int value) { - return ToStatus(XPRSsetintcontrol(xpress_model_, attribute, value)); -} - absl::StatusOr Xpress::GetDoubleAttr(int attribute) const { double result; RETURN_IF_ERROR(ToStatus(XPRSgetdblattrib(xpress_model_, attribute, &result))) @@ -191,7 +227,7 @@ absl::StatusOr> Xpress::GetPrimalValues() const { XPRSgetintattrib(xpress_model_, XPRS_COLS, &nVars); std::vector values(nVars); RETURN_IF_ERROR(ToStatus( - XPRSgetlpsol(xpress_model_, values.data(), nullptr, nullptr, nullptr))) + XPRSgetsolution(xpress_model_, nullptr, values.data(), 0, nVars - 1))) << "Error getting Xpress LP solution"; return values; } @@ -212,8 +248,8 @@ absl::StatusOr> Xpress::GetConstraintDuals() const { int nCons = GetNumberOfConstraints(); double values[nCons]; RETURN_IF_ERROR( - ToStatus(XPRSgetlpsol(xpress_model_, nullptr, nullptr, values, nullptr))) - << "Failed to retrieve LP solution from XPRESS"; + ToStatus(XPRSgetduals(xpress_model_, nullptr, values, 0, nCons - 1))) + << "Failed to retrieve duals from XPRESS"; std::vector result(values, values + nCons); return result; } @@ -221,7 +257,7 @@ absl::StatusOr> Xpress::GetReducedCostValues() const { int nVars = GetNumberOfVariables(); double values[nVars]; RETURN_IF_ERROR( - ToStatus(XPRSgetlpsol(xpress_model_, nullptr, nullptr, nullptr, values))) + ToStatus(XPRSgetredcosts(xpress_model_, nullptr, values, 0, nVars - 1))) << "Failed to retrieve LP solution from XPRESS"; std::vector result(values, values + nVars); return result; @@ -237,6 +273,10 @@ absl::Status Xpress::GetBasis(std::vector& rowBasis, absl::Status Xpress::SetStartingBasis(std::vector& rowBasis, std::vector& colBasis) const { + if (rowBasis.size() != colBasis.size()) { + return absl::InvalidArgumentError( + "Row basis and column basis must be of same size."); + } return ToStatus( XPRSloadbasis(xpress_model_, rowBasis.data(), colBasis.data())); } diff --git a/ortools/math_opt/solvers/xpress/g_xpress.h b/ortools/math_opt/solvers/xpress/g_xpress.h index a18ebabd1d2..9291db70430 100644 --- a/ortools/math_opt/solvers/xpress/g_xpress.h +++ b/ortools/math_opt/solvers/xpress/g_xpress.h @@ -30,6 +30,7 @@ #include #include #include +#include #include "absl/status/status.h" #include "absl/status/statusor.h" @@ -48,8 +49,11 @@ class Xpress { ~Xpress(); + absl::StatusOr GetIntControl(int control) const; + absl::Status SetIntControl(int control, int value); + absl::Status ResetIntControl(int control); // reset to default value + absl::StatusOr GetIntAttr(int attribute) const; - absl::Status SetIntAttr(int attribute, int value); absl::StatusOr GetDoubleAttr(int attribute) const; @@ -106,6 +110,9 @@ class Xpress { absl::Status ToStatus( int xprs_err, absl::StatusCode code = absl::StatusCode::kInvalidArgument) const; + + std::map int_control_defaults_; + void initIntControlDefaults(); }; } // namespace operations_research::math_opt diff --git a/ortools/math_opt/solvers/xpress_solver.cc b/ortools/math_opt/solvers/xpress_solver.cc index 7538841b542..a0d38907814 100644 --- a/ortools/math_opt/solvers/xpress_solver.cc +++ b/ortools/math_opt/solvers/xpress_solver.cc @@ -256,7 +256,7 @@ absl::Status XpressSolver::CallXpressSolve( const SolveParametersProto& parameters) { // Enable screen output right before solve if (parameters.enable_output()) { - RETURN_IF_ERROR(xpress_->SetIntAttr(XPRS_OUTPUTLOG, 1)) + RETURN_IF_ERROR(xpress_->SetIntControl(XPRS_OUTPUTLOG, 1)) << "Unable to enable XPRESS logs"; } // Solve @@ -275,7 +275,7 @@ absl::Status XpressSolver::CallXpressSolve( } // Disable screen output right after solve if (parameters.enable_output()) { - RETURN_IF_ERROR(xpress_->SetIntAttr(XPRS_OUTPUTLOG, 0)) + RETURN_IF_ERROR(xpress_->SetIntControl(XPRS_OUTPUTLOG, 0)) << "Unable to disable XPRESS logs"; } return absl::OkStatus(); @@ -285,17 +285,20 @@ absl::Status XpressSolver::SetLpIterLimits( const SolveParametersProto& parameters) { // If the user has set no limits, we still have to reset the limits // explicitly to their default values, else the parameters could be kept - // between solves. XPRESS default value for XPRS_LPITERLIMIT is 2147483647 - int lp_iter_limit = parameters.has_iteration_limit() - ? parameters.iteration_limit() - : 2147483647; - // XPRESS default value for XPRS_BARITERLIMIT is 500 - int bar_iter_limit = - parameters.has_iteration_limit() ? parameters.iteration_limit() : 500; - RETURN_IF_ERROR(xpress_->SetIntAttr(XPRS_LPITERLIMIT, lp_iter_limit)) - << "Could not set XPRS_LPITERLIMIT"; - RETURN_IF_ERROR(xpress_->SetIntAttr(XPRS_BARITERLIMIT, bar_iter_limit)) - << "Could not set XPRS_BARITERLIMIT"; + // between solves. + if (parameters.has_iteration_limit()) { + RETURN_IF_ERROR( + xpress_->SetIntControl(XPRS_LPITERLIMIT, parameters.iteration_limit())) + << "Could not set XPRS_LPITERLIMIT"; + RETURN_IF_ERROR( + xpress_->SetIntControl(XPRS_BARITERLIMIT, parameters.iteration_limit())) + << "Could not set XPRS_BARITERLIMIT"; + } else { + RETURN_IF_ERROR(xpress_->ResetIntControl(XPRS_LPITERLIMIT)) + << "Could not reset XPRS_LPITERLIMIT to its default value"; + RETURN_IF_ERROR(xpress_->ResetIntControl(XPRS_BARITERLIMIT)) + << "Could not reset XPRS_BARITERLIMIT to its default value"; + } return absl::OkStatus(); } @@ -371,7 +374,8 @@ absl::StatusOr XpressSolver::GetLpSolution( bool XpressSolver::isFeasible() const { if (is_mip_) { - return xpress_status_ == XPRS_MIP_OPTIMAL; + return xpress_status_ == XPRS_MIP_OPTIMAL || + xpress_status_ == XPRS_MIP_SOLUTION; } else { return xpress_status_ == XPRS_LP_OPTIMAL || xpress_status_ == XPRS_LP_UNFINISHED; diff --git a/ortools/xpress/environment.cc b/ortools/xpress/environment.cc index 13b9fccb4da..1bad91c1ff9 100644 --- a/ortools/xpress/environment.cc +++ b/ortools/xpress/environment.cc @@ -69,6 +69,9 @@ std::function XPRSgetrhsr std::function XPRSgetlb = nullptr; std::function XPRSgetub = nullptr; std::function XPRSgetcoef = nullptr; +std::function XPRSgetsolution = nullptr; +std::function XPRSgetduals = nullptr; +std::function XPRSgetredcosts = nullptr; std::function XPRSaddrows = nullptr; std::function XPRSdelrows = nullptr; std::function XPRSaddcols = nullptr; @@ -91,6 +94,7 @@ std::function XPRSgetmipsol = nu std::function XPRSchgobj = nullptr; std::function XPRSchgcoef = nullptr; std::function XPRSchgmcoef = nullptr; +std::function XPRSchgmcoef64 = nullptr; std::function XPRSchgrhs = nullptr; std::function XPRSchgrhsrange = nullptr; std::function XPRSchgrowtype = nullptr; @@ -99,6 +103,7 @@ std::function XPRSaddcbmessage = nullptr; std::function XPRSlpoptimize = nullptr; std::function XPRSmipoptimize = nullptr; +std::function XPRSoptimize = nullptr; void LoadXpressFunctions(DynamicLibrary* xpress_dynamic_library) { // This was generated with the parse_header_xpress.py script. @@ -133,6 +138,9 @@ void LoadXpressFunctions(DynamicLibrary* xpress_dynamic_library) { xpress_dynamic_library->GetFunction(&XPRSgetlb, "XPRSgetlb"); xpress_dynamic_library->GetFunction(&XPRSgetub, "XPRSgetub"); xpress_dynamic_library->GetFunction(&XPRSgetcoef, "XPRSgetcoef"); + xpress_dynamic_library->GetFunction(&XPRSgetsolution, "XPRSgetsolution"); + xpress_dynamic_library->GetFunction(&XPRSgetduals, "XPRSgetduals"); + xpress_dynamic_library->GetFunction(&XPRSgetredcosts, "XPRSgetredcosts"); xpress_dynamic_library->GetFunction(&XPRSaddrows, "XPRSaddrows"); xpress_dynamic_library->GetFunction(&XPRSdelrows, "XPRSdelrows"); xpress_dynamic_library->GetFunction(&XPRSaddcols, "XPRSaddcols"); @@ -155,6 +163,7 @@ void LoadXpressFunctions(DynamicLibrary* xpress_dynamic_library) { xpress_dynamic_library->GetFunction(&XPRSchgobj, "XPRSchgobj"); xpress_dynamic_library->GetFunction(&XPRSchgcoef, "XPRSchgcoef"); xpress_dynamic_library->GetFunction(&XPRSchgmcoef, "XPRSchgmcoef"); + xpress_dynamic_library->GetFunction(&XPRSchgmcoef64, "XPRSchgmcoef64"); xpress_dynamic_library->GetFunction(&XPRSchgrhs, "XPRSchgrhs"); xpress_dynamic_library->GetFunction(&XPRSchgrhsrange, "XPRSchgrhsrange"); xpress_dynamic_library->GetFunction(&XPRSchgrowtype, "XPRSchgrowtype"); @@ -163,6 +172,7 @@ void LoadXpressFunctions(DynamicLibrary* xpress_dynamic_library) { xpress_dynamic_library->GetFunction(&XPRSaddcbmessage, "XPRSaddcbmessage"); xpress_dynamic_library->GetFunction(&XPRSlpoptimize, "XPRSlpoptimize"); xpress_dynamic_library->GetFunction(&XPRSmipoptimize, "XPRSmipoptimize"); + xpress_dynamic_library->GetFunction(&XPRSoptimize, "XPRSoptimize"); } // clang-format on diff --git a/ortools/xpress/environment.h b/ortools/xpress/environment.h index ab5cc400d6a..a1015e769fd 100644 --- a/ortools/xpress/environment.h +++ b/ortools/xpress/environment.h @@ -492,6 +492,9 @@ extern std::function XPRS extern std::function XPRSgetlb; extern std::function XPRSgetub; extern std::function XPRSgetcoef; +extern std::function XPRSgetsolution; +extern std::function XPRSgetduals; +extern std::function XPRSgetredcosts; extern std::function XPRSaddrows; extern std::function XPRSdelrows; extern std::function XPRSaddcols; @@ -514,6 +517,7 @@ extern std::function XPRSgetmips extern std::function XPRSchgobj; extern std::function XPRSchgcoef; extern std::function XPRSchgmcoef; +extern std::function XPRSchgmcoef64; extern std::function XPRSchgrhs; extern std::function XPRSchgrhsrange; extern std::function XPRSchgrowtype; @@ -522,6 +526,8 @@ extern std::function XPRSaddcbmessage; extern std::function XPRSlpoptimize; extern std::function XPRSmipoptimize; +extern std::function XPRSoptimize; + // clang-format on } // namespace operations_research From c9419e7ba45b6882dcfa30d1b9e4c5d7eb113114 Mon Sep 17 00:00:00 2001 From: Peter Mitri Date: Fri, 3 Jan 2025 17:18:21 +0100 Subject: [PATCH 13/14] remove redundant DCHECK --- ortools/math_opt/solvers/xpress/g_xpress.cc | 1 - 1 file changed, 1 deletion(-) diff --git a/ortools/math_opt/solvers/xpress/g_xpress.cc b/ortools/math_opt/solvers/xpress/g_xpress.cc index 1488de895de..0d1d8dfe4bc 100644 --- a/ortools/math_opt/solvers/xpress/g_xpress.cc +++ b/ortools/math_opt/solvers/xpress/g_xpress.cc @@ -65,7 +65,6 @@ absl::StatusOr> Xpress::New( CHECK(correctlyLoaded); XPRSprob model; CHECK_EQ(kXpressOk, XPRScreateprob(&model)); - DCHECK(model != nullptr); // should not be NULL if status=0 CHECK_EQ(kXpressOk, XPRSaddcbmessage(model, printXpressMessage, nullptr, 0)); return absl::WrapUnique(new Xpress(model)); } From 0ae4346e12c388d28b5a1dcf948bd5d95db61bfc Mon Sep 17 00:00:00 2001 From: Peter Mitri Date: Mon, 6 Jan 2025 17:12:11 +0100 Subject: [PATCH 14/14] compute dual bound --- ortools/math_opt/solvers/xpress_solver.cc | 84 +++++++++++++---------- ortools/math_opt/solvers/xpress_solver.h | 21 ++++-- 2 files changed, 63 insertions(+), 42 deletions(-) diff --git a/ortools/math_opt/solvers/xpress_solver.cc b/ortools/math_opt/solvers/xpress_solver.cc index a0d38907814..3274fc48995 100644 --- a/ortools/math_opt/solvers/xpress_solver.cc +++ b/ortools/math_opt/solvers/xpress_solver.cc @@ -151,8 +151,8 @@ absl::Status XpressSolver::AddNewLinearConstraints( char sense = XPRS_EQUAL; double rhs = 0.0; double rng = 0.0; - const bool lb_is_xprs_neg_inf = lb <= XPRS_MINUSINFINITY; - const bool ub_is_xprs_pos_inf = ub >= XPRS_PLUSINFINITY; + const bool lb_is_xprs_neg_inf = lb <= kMinusInf; + const bool ub_is_xprs_pos_inf = ub >= kPlusInf; if (lb_is_xprs_neg_inf && !ub_is_xprs_pos_inf) { sense = XPRS_LESS_EQUAL; rhs = ub; @@ -232,8 +232,9 @@ absl::StatusOr XpressSolver::Solve( RETURN_IF_ERROR(CallXpressSolve(parameters)) << "Error during XPRESS solve"; - ASSIGN_OR_RETURN(SolveResultProto solve_result, - ExtractSolveResultProto(start, model_parameters)); + ASSIGN_OR_RETURN( + SolveResultProto solve_result, + ExtractSolveResultProto(start, model_parameters, parameters)); return solve_result; } @@ -303,16 +304,19 @@ absl::Status XpressSolver::SetLpIterLimits( } absl::StatusOr XpressSolver::ExtractSolveResultProto( - absl::Time start, const ModelSolveParametersProto& model_parameters) { + absl::Time start, const ModelSolveParametersProto& model_parameters, + const SolveParametersProto& solve_parameters) { SolveResultProto result; ASSIGN_OR_RETURN(SolutionsAndClaims solution_and_claims, - GetSolutions(model_parameters)); + GetSolutions(model_parameters, solve_parameters)); for (SolutionProto& solution : solution_and_claims.solutions) { *result.add_solutions() = std::move(solution); } ASSIGN_OR_RETURN(*result.mutable_solve_stats(), GetSolveStats(start)); - ASSIGN_OR_RETURN(const double best_primal_bound, GetBestPrimalBound()); - ASSIGN_OR_RETURN(const double best_dual_bound, GetBestDualBound()); + ASSIGN_OR_RETURN(const double best_primal_bound, + GetBestPrimalBound(solve_parameters)); + ASSIGN_OR_RETURN(const double best_dual_bound, + GetBestDualBound(solve_parameters)); auto solution_claims = solution_and_claims.solution_claims; ASSIGN_OR_RETURN(*result.mutable_termination(), ConvertTerminationReason(solution_claims, best_primal_bound, @@ -320,30 +324,46 @@ absl::StatusOr XpressSolver::ExtractSolveResultProto( return result; } -absl::StatusOr XpressSolver::GetBestPrimalBound() const { - return xpress_->GetDoubleAttr(is_mip_ ? XPRS_MIPOBJVAL : XPRS_LPOBJVAL); +absl::StatusOr XpressSolver::GetBestPrimalBound( + const SolveParametersProto& parameters) const { + // TODO: handle MIP + if (parameters.lp_algorithm() == LP_ALGORITHM_PRIMAL_SIMPLEX && + isFeasible() || + xpress_status_ == XPRS_LP_OPTIMAL) { + return xpress_->GetDoubleAttr(XPRS_LPOBJVAL); + } + return is_maximize_ ? kMinusInf : kPlusInf; } -absl::StatusOr XpressSolver::GetBestDualBound() const { - // TODO: setting LP primal value as best dual bound. Can this be improved? - return xpress_->GetDoubleAttr(XPRS_LPOBJVAL); +absl::StatusOr XpressSolver::GetBestDualBound( + const SolveParametersProto& parameters) const { + // TODO: handle MIP + if (parameters.lp_algorithm() == LP_ALGORITHM_DUAL_SIMPLEX && isFeasible() || + xpress_status_ == XPRS_LP_OPTIMAL) { + return xpress_->GetDoubleAttr(XPRS_LPOBJVAL); + } + return is_maximize_ ? kPlusInf : kMinusInf; } absl::StatusOr XpressSolver::GetSolutions( - const ModelSolveParametersProto& model_parameters) { + const ModelSolveParametersProto& model_parameters, + const SolveParametersProto& solve_parameters) { if (is_mip_) { return absl::UnimplementedError("XpressSolver does not handle MIPs yet"); } else { - return GetLpSolution(model_parameters); + return GetLpSolution(model_parameters, solve_parameters); } } absl::StatusOr XpressSolver::GetLpSolution( - const ModelSolveParametersProto& model_parameters) { - ASSIGN_OR_RETURN(auto primal_solution_and_claim, - GetConvexPrimalSolutionIfAvailable(model_parameters)); - ASSIGN_OR_RETURN(auto dual_solution_and_claim, - GetConvexDualSolutionIfAvailable(model_parameters)); + const ModelSolveParametersProto& model_parameters, + const SolveParametersProto& solve_parameters) { + ASSIGN_OR_RETURN( + auto primal_solution_and_claim, + GetConvexPrimalSolutionIfAvailable(model_parameters, solve_parameters)); + ASSIGN_OR_RETURN( + auto dual_solution_and_claim, + GetConvexDualSolutionIfAvailable(model_parameters, solve_parameters)); ASSIGN_OR_RETURN(auto basis, GetBasisIfAvailable()); const SolutionClaims solution_claims = { .primal_feasible_solution_exists = @@ -406,14 +426,15 @@ SolutionStatusProto XpressSolver::getLpSolutionStatus() const { absl::StatusOr> XpressSolver::GetConvexPrimalSolutionIfAvailable( - const ModelSolveParametersProto& model_parameters) const { + const ModelSolveParametersProto& model_parameters, + const SolveParametersProto& solve_parameters) const { if (!isFeasible()) { return SolutionAndClaim{ .solution = std::nullopt, .feasible_solution_exists = false}; } PrimalSolutionProto primal_solution; primal_solution.set_feasibility_status(getLpSolutionStatus()); - ASSIGN_OR_RETURN(const double sol_val, GetBestPrimalBound()); + ASSIGN_OR_RETURN(const double sol_val, GetBestPrimalBound(solve_parameters)); primal_solution.set_objective_value(sol_val); XpressVectorToSparseDoubleVector(xpress_->GetPrimalValues().value(), variables_map_, @@ -428,7 +449,8 @@ XpressSolver::GetConvexPrimalSolutionIfAvailable( absl::StatusOr> XpressSolver::GetConvexDualSolutionIfAvailable( - const ModelSolveParametersProto& model_parameters) const { + const ModelSolveParametersProto& model_parameters, + const SolveParametersProto& solve_parameters) const { if (!isFeasible()) { return SolutionAndClaim{ .solution = std::nullopt, .feasible_solution_exists = false}; @@ -447,20 +469,12 @@ XpressSolver::GetConvexDualSolutionIfAvailable( XpressVectorToSparseDoubleVector(xprs_reduced_cost_values, variables_map_, *dual_solution.mutable_reduced_costs(), model_parameters.reduced_costs_filter()); - ASSIGN_OR_RETURN(const double sol_val, GetBestPrimalBound()); + ASSIGN_OR_RETURN(const double sol_val, GetBestDualBound(solve_parameters)); dual_solution.set_objective_value(sol_val); dual_solution.set_feasibility_status(getLpSolutionStatus()); - bool dual_feasible_solution_exists = - (dual_solution.feasibility_status() == SOLUTION_STATUS_FEASIBLE); - ASSIGN_OR_RETURN(const double best_dual_bound, GetBestDualBound()); - // const double best_dual_bound = is_maximize_ ? kMinusInf : kPlusInf; - if (dual_feasible_solution_exists || std::isfinite(best_dual_bound)) { - dual_feasible_solution_exists = true; - } else if (xpress_status_ == XPRS_LP_OPTIMAL) { - return absl::InternalError( - "Xpress status is XPRS_LP_OPTIMAL, but XPRS_BESTBOUND is " - "unavailable or infinite, and no dual feasible solution is returned"); - } + ASSIGN_OR_RETURN(const double best_dual_bound, + GetBestDualBound(solve_parameters)); + bool dual_feasible_solution_exists = std::isfinite(best_dual_bound); return SolutionAndClaim{ .solution = dual_solution, .feasible_solution_exists = dual_feasible_solution_exists}; diff --git a/ortools/math_opt/solvers/xpress_solver.h b/ortools/math_opt/solvers/xpress_solver.h index 0ce9cf31c94..33e35d0cda1 100644 --- a/ortools/math_opt/solvers/xpress_solver.h +++ b/ortools/math_opt/solvers/xpress_solver.h @@ -129,13 +129,17 @@ class XpressSolver : public SolverInterface { }; absl::StatusOr ExtractSolveResultProto( - absl::Time start, const ModelSolveParametersProto& model_parameters); + absl::Time start, const ModelSolveParametersProto& model_parameters, + const SolveParametersProto& solve_parameters); absl::StatusOr GetSolutions( - const ModelSolveParametersProto& model_parameters); + const ModelSolveParametersProto& model_parameters, + const SolveParametersProto& solve_parameters); absl::StatusOr GetSolveStats(absl::Time start) const; - absl::StatusOr GetBestDualBound() const; - absl::StatusOr GetBestPrimalBound() const; + absl::StatusOr GetBestPrimalBound( + const SolveParametersProto& parameters) const; + absl::StatusOr GetBestDualBound( + const SolveParametersProto& parameters) const; absl::StatusOr ConvertTerminationReason( SolutionClaims solution_claims, double best_primal_bound, @@ -144,16 +148,19 @@ class XpressSolver : public SolverInterface { // Returns solution information appropriate and available for an LP (linear // constraints + linear objective, only). absl::StatusOr GetLpSolution( - const ModelSolveParametersProto& model_parameters); + const ModelSolveParametersProto& model_parameters, + const SolveParametersProto& solve_parameters); bool isFeasible() const; // return bool field should be true if a primal solution exists. absl::StatusOr> GetConvexPrimalSolutionIfAvailable( - const ModelSolveParametersProto& model_parameters) const; + const ModelSolveParametersProto& model_parameters, + const SolveParametersProto& solve_parameters) const; absl::StatusOr> GetConvexDualSolutionIfAvailable( - const ModelSolveParametersProto& model_parameters) const; + const ModelSolveParametersProto& model_parameters, + const SolveParametersProto& solve_parameters) const; absl::StatusOr> GetBasisIfAvailable(); absl::Status AddNewLinearConstraints(