From 4d7c19cd9bf03d44ac6e39e274f028b53534a896 Mon Sep 17 00:00:00 2001 From: Corentin Le Molgat Date: Mon, 19 Sep 2022 18:30:50 +0200 Subject: [PATCH 1/2] pdlp: update --- ortools/pdlp/BUILD.bazel | 1 + ortools/pdlp/iteration_stats.cc | 7 +- ortools/pdlp/iteration_stats.h | 8 +- ortools/pdlp/iteration_stats_test.cc | 8 +- ortools/pdlp/primal_dual_hybrid_gradient.cc | 26 +- ortools/pdlp/primal_dual_hybrid_gradient.h | 2 - .../pdlp/primal_dual_hybrid_gradient_test.cc | 98 +++-- ortools/pdlp/quadratic_program.cc | 3 +- ortools/pdlp/quadratic_program.h | 1 - ortools/pdlp/quadratic_program_io.cc | 1 + ortools/pdlp/quadratic_program_test.cc | 6 +- ortools/pdlp/sharded_optimization_utils.cc | 2 +- .../pdlp/sharded_optimization_utils_test.cc | 9 +- ortools/pdlp/sharded_quadratic_program.cc | 1 + ortools/pdlp/sharder.cc | 1 + ortools/pdlp/sharder.h | 2 +- ortools/pdlp/solvers.proto | 76 +++- ortools/pdlp/solvers_proto_validation.cc | 109 ++++-- ortools/pdlp/solvers_proto_validation_test.cc | 146 +++++-- ortools/pdlp/termination.cc | 96 ++++- ortools/pdlp/termination.h | 11 +- ortools/pdlp/termination_test.cc | 355 ++++++++++++++++-- ortools/pdlp/test_util.cc | 1 + ortools/pdlp/test_util.h | 4 +- ortools/pdlp/test_util_test.cc | 2 +- ortools/pdlp/trust_region.cc | 3 +- ortools/pdlp/trust_region.h | 2 +- ortools/pdlp/trust_region_test.cc | 3 + 28 files changed, 792 insertions(+), 192 deletions(-) diff --git a/ortools/pdlp/BUILD.bazel b/ortools/pdlp/BUILD.bazel index 0edb56a476c..d581bdb3d08 100644 --- a/ortools/pdlp/BUILD.bazel +++ b/ortools/pdlp/BUILD.bazel @@ -267,6 +267,7 @@ cc_test( ":gtest_main", ":solvers_cc_proto", ":solvers_proto_validation", + "//ortools/base:protobuf_util", "@com_google_absl//absl/status", ], ) diff --git a/ortools/pdlp/iteration_stats.cc b/ortools/pdlp/iteration_stats.cc index 29635e4a1cd..359a0760bd6 100644 --- a/ortools/pdlp/iteration_stats.cc +++ b/ortools/pdlp/iteration_stats.cc @@ -16,20 +16,15 @@ #include #include #include -#include #include #include -#include #include #include #include "Eigen/Core" #include "Eigen/SparseCore" #include "absl/random/distributions.h" -#include "absl/strings/str_cat.h" -#include "absl/strings/str_format.h" -#include "absl/strings/string_view.h" -#include "ortools/base/logging.h" +#include "ortools/base/check.h" #include "ortools/base/mathutil.h" #include "ortools/pdlp/quadratic_program.h" #include "ortools/pdlp/sharded_quadratic_program.h" diff --git a/ortools/pdlp/iteration_stats.h b/ortools/pdlp/iteration_stats.h index 6382dcc9bb5..b3bce16b9d6 100644 --- a/ortools/pdlp/iteration_stats.h +++ b/ortools/pdlp/iteration_stats.h @@ -14,9 +14,7 @@ #ifndef PDLP_ITERATION_STATS_H_ #define PDLP_ITERATION_STATS_H_ -#include #include -#include #include #include "Eigen/Core" @@ -79,17 +77,17 @@ Eigen::VectorXd ReducedCosts(const ShardedQuadraticProgram& scaled_sharded_qp, bool use_zero_primal_objective = false); // Finds and returns the ConvergenceInformation with the specified -// candidate_type, or absl::nullopt if no such candidate exists. +// candidate_type, or std::nullopt if no such candidate exists. std::optional GetConvergenceInformation( const IterationStats& stats, PointType candidate_type); // Finds and returns the InfeasibilityInformation with the specified -// candidate_type, or absl::nullopt if no such candidate exists. +// candidate_type, or std::nullopt if no such candidate exists. std::optional GetInfeasibilityInformation( const IterationStats& stats, PointType candidate_type); // Finds and returns the PointMetadata with the specified -// point_type, or absl::nullopt if no such point exists. +// point_type, or std::nullopt if no such point exists. std::optional GetPointMetadata(const IterationStats& stats, PointType point_type); diff --git a/ortools/pdlp/iteration_stats_test.cc b/ortools/pdlp/iteration_stats_test.cc index f50f98a2e1b..1e5236dbb1a 100644 --- a/ortools/pdlp/iteration_stats_test.cc +++ b/ortools/pdlp/iteration_stats_test.cc @@ -14,7 +14,7 @@ #include "ortools/pdlp/iteration_stats.h" #include -#include +#include #include #include "Eigen/Core" @@ -218,7 +218,7 @@ TEST(GetConvergenceInformation, GetsCorrectEntry) { EXPECT_THAT( GetConvergenceInformation(test_stats, POINT_TYPE_ITERATE_DIFFERENCE), - Eq(absl::nullopt)); + Eq(std::nullopt)); } TEST(GetInfeasibilityInformation, GetsCorrectEntry) { @@ -246,7 +246,7 @@ TEST(GetInfeasibilityInformation, GetsCorrectEntry) { EXPECT_THAT( GetInfeasibilityInformation(test_stats, POINT_TYPE_ITERATE_DIFFERENCE), - Eq(absl::nullopt)); + Eq(std::nullopt)); } TEST(GetPointMetadata, GetsCorrectEntry) { @@ -273,7 +273,7 @@ TEST(GetPointMetadata, GetsCorrectEntry) { EXPECT_EQ(current_info->active_primal_variable_count(), 1); EXPECT_THAT(GetPointMetadata(test_stats, POINT_TYPE_ITERATE_DIFFERENCE), - Eq(absl::nullopt)); + Eq(std::nullopt)); } } // namespace diff --git a/ortools/pdlp/primal_dual_hybrid_gradient.cc b/ortools/pdlp/primal_dual_hybrid_gradient.cc index 501d6ff4daf..080e1bad6c1 100644 --- a/ortools/pdlp/primal_dual_hybrid_gradient.cc +++ b/ortools/pdlp/primal_dual_hybrid_gradient.cc @@ -174,9 +174,9 @@ std::string ToString(const IterationStats& iter_stats, } if (convergence_information.has_value()) { const RelativeConvergenceInformation relative_information = - ComputeRelativeResiduals(termination_criteria.eps_optimal_absolute(), - termination_criteria.eps_optimal_relative(), - bound_norms, *convergence_information); + ComputeRelativeResiduals( + EffectiveOptimalityCriteria(termination_criteria), bound_norms, + *convergence_information); return absl::StrCat(iteration_string, " | ", ToString(*convergence_information, relative_information, termination_criteria.optimality_norm())); @@ -199,9 +199,9 @@ std::string ToShortString(const IterationStats& iter_stats, } if (convergence_information.has_value()) { const RelativeConvergenceInformation relative_information = - ComputeRelativeResiduals(termination_criteria.eps_optimal_absolute(), - termination_criteria.eps_optimal_relative(), - bound_norms, *convergence_information); + ComputeRelativeResiduals( + EffectiveOptimalityCriteria(termination_criteria), bound_norms, + *convergence_information); return absl::StrCat( iteration_string, " | ", ToShortString(*convergence_information, relative_information, @@ -1273,7 +1273,7 @@ std::optional Solver::MajorIterationAndTerminationCheck( *solve_log.add_iteration_stats() = stats; } ApplyRestartChoice(restart); - return absl::nullopt; + return std::nullopt; } void Solver::ResetAverageToCurrent() { @@ -1589,19 +1589,19 @@ std::optional Solver::ApplyPresolveIfEnabled( std::optional* const initial_solution) { const bool presolve_enabled = params_.presolve_options().use_glop(); if (!presolve_enabled) { - return absl::nullopt; + return std::nullopt; } if (!IsLinearProgram(WorkingQp())) { LOG(WARNING) << "Skipping presolve, which is only supported for linear programs"; - return absl::nullopt; + return std::nullopt; } absl::StatusOr model = QpToMpModelProto(WorkingQp()); if (!model.ok()) { LOG(WARNING) << "Skipping presolve because of error converting to MPModelProto: " << model.status(); - return absl::nullopt; + return std::nullopt; } if (initial_solution->has_value()) { LOG(WARNING) << "Ignoring initial solution. Initial solutions " @@ -1643,7 +1643,7 @@ std::optional Solver::ApplyPresolveIfEnabled( row_scaling_vec_ = OnesVector(sharded_working_qp_.DualSharder()); return GlopStatusToTerminationReason(presolve_info_->preprocessor.status()); } - return absl::nullopt; + return std::nullopt; } PrimalAndDualSolution Solver::RecoverOriginalSolution( @@ -1933,7 +1933,7 @@ SolverResult Solver::Solve( double inverse_step_size; const auto lipschitz_result = EstimateMaximumSingularValueOfConstraintMatrix( - sharded_working_qp_, absl::nullopt, absl::nullopt, + sharded_working_qp_, std::nullopt, std::nullopt, /*desired_relative_error=*/0.2, /*failure_probability=*/0.0005, random); // With high probability, the estimate of the lipschitz term is within @@ -2032,7 +2032,7 @@ SolverResult PrimalDualHybridGradient( QuadraticProgram qp, const PrimalDualHybridGradientParams& params, const std::atomic* interrupt_solve, IterationStatsCallback iteration_stats_callback) { - return PrimalDualHybridGradient(std::move(qp), params, absl::nullopt, + return PrimalDualHybridGradient(std::move(qp), params, std::nullopt, interrupt_solve, std::move(iteration_stats_callback)); } diff --git a/ortools/pdlp/primal_dual_hybrid_gradient.h b/ortools/pdlp/primal_dual_hybrid_gradient.h index 5f50054dc0c..b4e2606dfba 100644 --- a/ortools/pdlp/primal_dual_hybrid_gradient.h +++ b/ortools/pdlp/primal_dual_hybrid_gradient.h @@ -19,8 +19,6 @@ #include #include "Eigen/Core" -#include "absl/status/statusor.h" -#include "ortools/linear_solver/linear_solver.pb.h" #include "ortools/lp_data/lp_data.h" #include "ortools/pdlp/quadratic_program.h" #include "ortools/pdlp/solve_log.pb.h" diff --git a/ortools/pdlp/primal_dual_hybrid_gradient_test.cc b/ortools/pdlp/primal_dual_hybrid_gradient_test.cc index 02a87deaca1..60ac4ed752d 100644 --- a/ortools/pdlp/primal_dual_hybrid_gradient_test.cc +++ b/ortools/pdlp/primal_dual_hybrid_gradient_test.cc @@ -18,6 +18,7 @@ #include #include #include +#include #include #include #include @@ -30,6 +31,7 @@ #include "absl/strings/str_cat.h" #include "gmock/gmock.h" #include "gtest/gtest.h" +#include "ortools/base/check.h" #include "ortools/base/logging.h" #include "ortools/glop/parameters.pb.h" #include "ortools/linear_solver/linear_solver.pb.h" @@ -75,14 +77,19 @@ PrimalDualHybridGradientParams CreateSolverParams( PrimalDualHybridGradientParams::MALITSKY_POCK_LINESEARCH_RULE); } - params.mutable_termination_criteria()->set_eps_optimal_relative(0.0); + params.mutable_termination_criteria() + ->mutable_simple_optimality_criteria() + ->set_eps_optimal_relative(0.0); if (use_iteration_limit) { // This effectively forces convergence on the iteration limit only. params.mutable_termination_criteria()->set_iteration_limit(iteration_limit); - params.mutable_termination_criteria()->set_eps_optimal_absolute(0.0); + params.mutable_termination_criteria() + ->mutable_simple_optimality_criteria() + ->set_eps_optimal_absolute(0.0); } else { - params.mutable_termination_criteria()->set_eps_optimal_absolute( - eps_optimal_absolute); + params.mutable_termination_criteria() + ->mutable_simple_optimality_criteria() + ->set_eps_optimal_absolute(eps_optimal_absolute); } if (use_diagonal_qp_trust_region_solver) { params.set_use_diagonal_qp_trust_region_solver(true); @@ -613,8 +620,12 @@ PrimalDualHybridGradientParams ParamsWithNoLimits() { PrimalDualHybridGradientParams params; // This disables the termination limits. A termination criteria must be set // for the solver to terminate. - params.mutable_termination_criteria()->set_eps_optimal_relative(0.0); - params.mutable_termination_criteria()->set_eps_optimal_absolute(0.0); + params.mutable_termination_criteria() + ->mutable_simple_optimality_criteria() + ->set_eps_optimal_relative(0.0); + params.mutable_termination_criteria() + ->mutable_simple_optimality_criteria() + ->set_eps_optimal_absolute(0.0); params.set_record_iteration_stats(true); return params; } @@ -882,24 +893,24 @@ TEST(PrimalDualHybridGradientTest, EXPECT_EQ(output.solve_log.iteration_stats().size(), iteration_limit + 1); for (const auto& stats : output.solve_log.iteration_stats()) { EXPECT_NE(GetConvergenceInformation(stats, POINT_TYPE_CURRENT_ITERATE), - absl::nullopt); + std::nullopt); EXPECT_NE(GetInfeasibilityInformation(stats, POINT_TYPE_CURRENT_ITERATE), - absl::nullopt); + std::nullopt); EXPECT_NE(GetPointMetadata(stats, POINT_TYPE_CURRENT_ITERATE), - absl::nullopt); + std::nullopt); if (stats.iteration_number() > 0) { EXPECT_NE(GetConvergenceInformation(stats, POINT_TYPE_AVERAGE_ITERATE), - absl::nullopt); + std::nullopt); EXPECT_NE(GetInfeasibilityInformation(stats, POINT_TYPE_AVERAGE_ITERATE), - absl::nullopt); + std::nullopt); EXPECT_NE(GetPointMetadata(stats, POINT_TYPE_AVERAGE_ITERATE), - absl::nullopt); + std::nullopt); EXPECT_NE( GetInfeasibilityInformation(stats, POINT_TYPE_ITERATE_DIFFERENCE), - absl::nullopt); + std::nullopt); EXPECT_NE(GetPointMetadata(stats, POINT_TYPE_ITERATE_DIFFERENCE), - absl::nullopt); + std::nullopt); } } } @@ -1214,7 +1225,9 @@ PrimalAndDualSolution TinyLpSolution() { TEST(PrimalDualHybridGradientTest, WarmStartedAtOptimum) { PrimalDualHybridGradientParams params = ParamsWithNoLimits(); - params.mutable_termination_criteria()->set_eps_optimal_absolute(1.0e-10); + params.mutable_termination_criteria() + ->mutable_simple_optimality_criteria() + ->set_eps_optimal_absolute(1.0e-10); SolverResult output = PrimalDualHybridGradient(TinyLp(), params, TinyLpSolution()); @@ -1260,8 +1273,12 @@ TEST(PrimalDualHybridGradientTest, EmptyQp) { TEST(PrimalDualHybridGradientTest, RespectsInterrupt) { std::atomic interrupt_solve; PrimalDualHybridGradientParams params; - params.mutable_termination_criteria()->set_eps_optimal_absolute(0.0); - params.mutable_termination_criteria()->set_eps_optimal_relative(0.0); + params.mutable_termination_criteria() + ->mutable_simple_optimality_criteria() + ->set_eps_optimal_absolute(0.0); + params.mutable_termination_criteria() + ->mutable_simple_optimality_criteria() + ->set_eps_optimal_relative(0.0); interrupt_solve.store(true); const SolverResult output = @@ -1273,8 +1290,12 @@ TEST(PrimalDualHybridGradientTest, RespectsInterrupt) { TEST(PrimalDualHybridGradientTest, RespectsInterruptFromCallback) { std::atomic interrupt_solve; PrimalDualHybridGradientParams params; - params.mutable_termination_criteria()->set_eps_optimal_absolute(0.0); - params.mutable_termination_criteria()->set_eps_optimal_relative(0.0); + params.mutable_termination_criteria() + ->mutable_simple_optimality_criteria() + ->set_eps_optimal_absolute(0.0); + params.mutable_termination_criteria() + ->mutable_simple_optimality_criteria() + ->set_eps_optimal_relative(0.0); interrupt_solve.store(false); auto callback = [&](const IterationCallbackInfo& info) { @@ -1293,8 +1314,12 @@ TEST(PrimalDualHybridGradientTest, RespectsInterruptFromCallback) { TEST(PrimalDualHybridGradientTest, IgnoresFalseInterrupt) { std::atomic interrupt_solve; PrimalDualHybridGradientParams params; - params.mutable_termination_criteria()->set_eps_optimal_absolute(0.0); - params.mutable_termination_criteria()->set_eps_optimal_relative(0.0); + params.mutable_termination_criteria() + ->mutable_simple_optimality_criteria() + ->set_eps_optimal_absolute(0.0); + params.mutable_termination_criteria() + ->mutable_simple_optimality_criteria() + ->set_eps_optimal_relative(0.0); params.mutable_termination_criteria()->set_kkt_matrix_pass_limit(1); interrupt_solve.store(false); @@ -1328,6 +1353,37 @@ TEST(PrimalDualHybridGradientTest, HugeNumShards) { TERMINATION_REASON_ITERATION_LIMIT); } +TEST(PrimalDualHybridGradientTest, DetailedTerminationCriteria) { + const int iteration_upperbound = 300; + PrimalDualHybridGradientParams params = CreateSolverParams( + iteration_upperbound, + /*eps_optimal_absolute=*/1.0e-5, /*enable_scaling=*/true, + /*num_threads=*/4, /*use_iteration_limit=*/false, + /*use_malitsky_pock_linesearch=*/false, + /*use_diagonal_qp_trust_region_solver=*/false); + params.set_major_iteration_frequency(60); + params.mutable_termination_criteria()->clear_simple_optimality_criteria(); + auto* opt_criteria = params.mutable_termination_criteria() + ->mutable_detailed_optimality_criteria(); + opt_criteria->set_eps_optimal_primal_residual_absolute(1e-5); + opt_criteria->set_eps_optimal_primal_residual_relative(0.0); + opt_criteria->set_eps_optimal_dual_residual_absolute(1e-5); + opt_criteria->set_eps_optimal_dual_residual_relative(0.0); + opt_criteria->set_eps_optimal_objective_gap_absolute(1e-5); + opt_criteria->set_eps_optimal_objective_gap_relative(0.0); + + SolverResult output = PrimalDualHybridGradient(TinyLp(), params); + VerifyTerminationReasonAndIterationCount(params, output, + /*use_iteration_limit=*/false); + VerifyObjectiveValues(output, -1.0, 1.0e-4); + EXPECT_THAT(output.primal_solution, + EigenArrayNear({1, 0, 6, 2}, 1.0e-4)); + EXPECT_THAT(output.dual_solution, + EigenArrayNear({0.5, 4.0, 0.0}, 1.0e-4)); + EXPECT_THAT(output.reduced_costs, + EigenArrayNear({0.0, 1.5, -3.5, 0.0}, 1.0e-4)); +} + // Verifies that the primal and dual solution satisfy the bounds constraints. // This function uses ASSERTS rather than EXPECTs because it's used with large // QPs that would spam the logs otherwise. diff --git a/ortools/pdlp/quadratic_program.cc b/ortools/pdlp/quadratic_program.cc index 1c1e503de74..6ff3642aa39 100644 --- a/ortools/pdlp/quadratic_program.cc +++ b/ortools/pdlp/quadratic_program.cc @@ -15,7 +15,6 @@ #include #include -#include #include #include #include @@ -27,7 +26,7 @@ #include "absl/status/status.h" #include "absl/status/statusor.h" #include "absl/strings/str_cat.h" -#include "ortools/base/logging.h" +#include "ortools/base/check.h" #include "ortools/base/status_macros.h" #include "ortools/linear_solver/linear_solver.pb.h" diff --git a/ortools/pdlp/quadratic_program.h b/ortools/pdlp/quadratic_program.h index d74bc352835..370da6197c9 100644 --- a/ortools/pdlp/quadratic_program.h +++ b/ortools/pdlp/quadratic_program.h @@ -25,7 +25,6 @@ #include "Eigen/SparseCore" #include "absl/status/status.h" #include "absl/status/statusor.h" -#include "ortools/base/logging.h" #include "ortools/linear_solver/linear_solver.pb.h" namespace operations_research::pdlp { diff --git a/ortools/pdlp/quadratic_program_io.cc b/ortools/pdlp/quadratic_program_io.cc index 797de089895..cab6ec50576 100644 --- a/ortools/pdlp/quadratic_program_io.cc +++ b/ortools/pdlp/quadratic_program_io.cc @@ -29,6 +29,7 @@ #include "absl/status/statusor.h" #include "absl/strings/match.h" #include "ortools/base/basictypes.h" +#include "ortools/base/check.h" #include "ortools/base/file.h" #include "ortools/base/helpers.h" #include "ortools/base/logging.h" diff --git a/ortools/pdlp/quadratic_program_test.cc b/ortools/pdlp/quadratic_program_test.cc index 62c244de7cb..00c4a7ca1a8 100644 --- a/ortools/pdlp/quadratic_program_test.cc +++ b/ortools/pdlp/quadratic_program_test.cc @@ -389,9 +389,9 @@ TEST(QpFromMpModelProtoTest, DoesNotIncludeNames) { QpFromMpModelProto(TinyModelWithNames(), /*relax_integer_variables=*/true, /*include_names=*/false); ASSERT_TRUE(lp.ok()) << lp.status(); - EXPECT_EQ(lp->problem_name, absl::nullopt); - EXPECT_EQ(lp->variable_names, absl::nullopt); - EXPECT_EQ(lp->constraint_names, absl::nullopt); + EXPECT_EQ(lp->problem_name, std::nullopt); + EXPECT_EQ(lp->variable_names, std::nullopt); + EXPECT_EQ(lp->constraint_names, std::nullopt); } TEST(QpFromMpModelProtoTest, IncludesNames) { diff --git a/ortools/pdlp/sharded_optimization_utils.cc b/ortools/pdlp/sharded_optimization_utils.cc index 29b657f0823..dcb9c3fa91f 100644 --- a/ortools/pdlp/sharded_optimization_utils.cc +++ b/ortools/pdlp/sharded_optimization_utils.cc @@ -16,7 +16,6 @@ #include #include #include -#include #include #include #include @@ -26,6 +25,7 @@ #include "Eigen/Core" #include "Eigen/SparseCore" #include "absl/random/distributions.h" +#include "ortools/base/check.h" #include "ortools/base/logging.h" #include "ortools/base/mathutil.h" #include "ortools/pdlp/quadratic_program.h" diff --git a/ortools/pdlp/sharded_optimization_utils_test.cc b/ortools/pdlp/sharded_optimization_utils_test.cc index 047d78e00dc..4dc40838c0f 100644 --- a/ortools/pdlp/sharded_optimization_utils_test.cc +++ b/ortools/pdlp/sharded_optimization_utils_test.cc @@ -15,6 +15,7 @@ #include #include +#include #include #include #include @@ -567,7 +568,7 @@ TEST(EstimateSingularValuesTest, CorrectForTestLp) { // The test_lp matrix is [ 2 1 1 2; 1 0 1 0; 4 0 0 0; 0 0 1.5 -1]. std::mt19937 random(1); auto result = EstimateMaximumSingularValueOfConstraintMatrix( - lp, absl::nullopt, absl::nullopt, + lp, std::nullopt, std::nullopt, /*desired_relative_error=*/0.01, /*failure_probability=*/0.001, random); EXPECT_NEAR(result.singular_value, 4.76945, 0.01); @@ -584,7 +585,7 @@ TEST(EstimateSingularValuesTest, CorrectForTestLpWithActivePrimalSubspace) { // so the projected matrix is [ 2 1 2; 1 1 0; 4 0 0; 0 1.5 -1]. std::mt19937 random(1); auto result = EstimateMaximumSingularValueOfConstraintMatrix( - lp, primal_solution, absl::nullopt, /*desired_relative_error=*/0.01, + lp, primal_solution, std::nullopt, /*desired_relative_error=*/0.01, /*failure_probability=*/0.001, random); EXPECT_NEAR(result.singular_value, 4.73818, 0.01); EXPECT_LT(result.num_iterations, 300); @@ -601,7 +602,7 @@ TEST(EstimateSingularValuesTest, CorrectForTestLpWithActiveDualSubspace) { // so the projected matrix is [ 2 1 1 2; 4 0 0 0; 0 0 1.5 -1]. std::mt19937 random(1); auto result = EstimateMaximumSingularValueOfConstraintMatrix( - lp, absl::nullopt, dual_solution, /*desired_relative_error=*/0.01, + lp, std::nullopt, dual_solution, /*desired_relative_error=*/0.01, /*failure_probability=*/0.001, random); EXPECT_NEAR(result.singular_value, 4.64203, 0.01); EXPECT_LT(result.num_iterations, 300); @@ -637,7 +638,7 @@ TEST(EstimateSingularValuesTest, CorrectForDiagonalLp) { // The diagonal_lp matrix is [ 2 0 0 0; 0 1 0 0; 0 0 -3 0; 0 0 0 -1]. std::mt19937 random(1); auto result = EstimateMaximumSingularValueOfConstraintMatrix( - lp, absl::nullopt, absl::nullopt, + lp, std::nullopt, std::nullopt, /*desired_relative_error=*/0.01, /*failure_probability=*/0.001, random); EXPECT_NEAR(result.singular_value, 3, 0.0001); diff --git a/ortools/pdlp/sharded_quadratic_program.cc b/ortools/pdlp/sharded_quadratic_program.cc index 5657e6fcf62..5159d551aa9 100644 --- a/ortools/pdlp/sharded_quadratic_program.cc +++ b/ortools/pdlp/sharded_quadratic_program.cc @@ -21,6 +21,7 @@ #include "Eigen/SparseCore" #include "absl/memory/memory.h" #include "absl/strings/string_view.h" +#include "ortools/base/check.h" #include "ortools/base/logging.h" #include "ortools/base/threadpool.h" #include "ortools/pdlp/quadratic_program.h" diff --git a/ortools/pdlp/sharder.cc b/ortools/pdlp/sharder.cc index 4212b0759e9..c457305a055 100644 --- a/ortools/pdlp/sharder.cc +++ b/ortools/pdlp/sharder.cc @@ -23,6 +23,7 @@ #include "Eigen/SparseCore" #include "absl/synchronization/blocking_counter.h" #include "absl/time/time.h" +#include "ortools/base/check.h" #include "ortools/base/logging.h" #include "ortools/base/mathutil.h" #include "ortools/base/threadpool.h" diff --git a/ortools/pdlp/sharder.h b/ortools/pdlp/sharder.h index 32ac79a83b4..d953cc4f5fa 100644 --- a/ortools/pdlp/sharder.h +++ b/ortools/pdlp/sharder.h @@ -21,7 +21,7 @@ #include "Eigen/Core" #include "Eigen/SparseCore" -#include "ortools/base/logging.h" +#include "ortools/base/check.h" #include "ortools/base/threadpool.h" namespace operations_research::pdlp { diff --git a/ortools/pdlp/solvers.proto b/ortools/pdlp/solvers.proto index 1afcfc4dd77..e220a6c9876 100644 --- a/ortools/pdlp/solvers.proto +++ b/ortools/pdlp/solvers.proto @@ -54,24 +54,70 @@ message TerminationCriteria { // If the algorithm terminates with termination_reason = // TERMINATION_REASON_OPTIMAL then // - // | primal_objective - dual_objective | <= eps_optimal_absolute + - // eps_optimal_relative * ( | primal_objective | + | dual_objective | ) - // norm(primal_residual, p) <= eps_optimal_absolute + eps_optimal_relative * - // norm(b, p) - // norm(dual_residual, p) <= eps_optimal_absolute + eps_optimal_relative * - // norm(c, p) + // When using DetailedOptimalityCriteria the conditions to declare a solution + // optimal are: + // | primal_objective - dual_objective | <= + // eps_optimal_objective_gap_absolute + + // eps_optimal_objective_gap_relative * + // ( | primal_objective | + | dual_objective | ) + // norm(primal_residual, p) <= + // eps_optimal_primal_residual_absolute + + // eps_optimal_primal_residual_relative * norm(b, p) + // norm(dual_residual, p) <= + // eps_optimal_dual_residual_absolute + + // eps_optimal_dual_residual_relative * norm(c, p) // It is possible to prove that a solution satisfying the above conditions - // also satisfies SCS's optimality conditions (see link above) with ϵ_pri = - // ϵ_dual = ϵ_gap = eps_optimal_absolute = eps_optimal_relative. (ϵ_pri, - // ϵ_dual, and ϵ_gap are SCS's parameters). + // also satisfies SCS's optimality conditions (see link above) with + // ϵ_pri = ϵ_dual = ϵ_gap = eps_optimal_*_absolute = eps_optimal_*_relative. + // (ϵ_pri, ϵ_dual, and ϵ_gap are SCS's parameters). + // When using SimpleOptimalityCriteria all the eps_optimal_*_absolute have the + // same value eps_optimal_absolute and all the eps_optimal_*_relative have the + // same value eps_optimal_relative. + + message SimpleOptimalityCriteria { + // Absolute tolerance on the primal residual, dual residual, and objective + // gap. + optional double eps_optimal_absolute = 1 [default = 1.0e-6]; + // Relative tolerance on the primal residual, dual residual, and objective + // gap. + optional double eps_optimal_relative = 2 [default = 1.0e-6]; + } + + message DetailedOptimalityCriteria { + // Absolute tolerance on the primal residual. + optional double eps_optimal_primal_residual_absolute = 1 [default = 1.0e-6]; + // Relative tolerance on the primal residual. + optional double eps_optimal_primal_residual_relative = 2 [default = 1.0e-6]; + // Absolute tolerance on the dual residual. + optional double eps_optimal_dual_residual_absolute = 3 [default = 1.0e-6]; + // Relative tolerance on the dual residual. + optional double eps_optimal_dual_residual_relative = 4 [default = 1.0e-6]; + // Absolute tolerance on the objective gap. + optional double eps_optimal_objective_gap_absolute = 5 [default = 1.0e-6]; + // Relative tolerance on the objective gap. + optional double eps_optimal_objective_gap_relative = 6 [default = 1.0e-6]; + } - // Absolute tolerance on the duality gap, primal feasibility, and dual - // feasibility. - optional double eps_optimal_absolute = 2 [default = 1.0e-6]; + // If none are set explicitly the deprecated eps_optimal_absolute and + // eps_optimal_relative are used to construct a SimpleOptimalityCriteria. + oneof optimality_criteria { + SimpleOptimalityCriteria simple_optimality_criteria = 9; + DetailedOptimalityCriteria detailed_optimality_criteria = 10; + } - // Relative tolerance on the duality gap, primal feasibility, and dual - // feasibility. - optional double eps_optimal_relative = 3 [default = 1.0e-6]; + // Absolute tolerance on primal residual, dual residual, and the objective + // gap. + // Deprecated, use simple_optimality_criteria instead. + // TODO(b/241462829) delete this deprecated field. + optional double eps_optimal_absolute = 2 + [deprecated = true, default = 1.0e-6]; + + // Relative tolerance on primal residual, dual residual, and the objective + // gap. + // Deprecated, use simple_optimality_criteria instead. + // TODO(b/241462829) delete this deprecated field. + optional double eps_optimal_relative = 3 + [deprecated = true, default = 1.0e-6]; // If the following two conditions hold we say that we have obtained an // approximate dual ray, which is an approximate certificate of primal diff --git a/ortools/pdlp/solvers_proto_validation.cc b/ortools/pdlp/solvers_proto_validation.cc index aeafc6eb12f..2a3ef693e8c 100644 --- a/ortools/pdlp/solvers_proto_validation.cc +++ b/ortools/pdlp/solvers_proto_validation.cc @@ -16,6 +16,8 @@ #include #include "absl/status/status.h" +#include "absl/strings/str_cat.h" +#include "absl/strings/string_view.h" #include "ortools/base/status_macros.h" #include "ortools/pdlp/solvers.pb.h" @@ -24,50 +26,85 @@ namespace operations_research::pdlp { using ::absl::InvalidArgumentError; using ::absl::OkStatus; +absl::Status CheckNonNegative(const double value, + const absl::string_view name) { + if (std::isnan(value)) { + return InvalidArgumentError(absl::StrCat(name, " is NAN")); + } + if (value < 0) { + return InvalidArgumentError(absl::StrCat(name, " must be non-negative")); + } + return absl::OkStatus(); +} absl::Status ValidateTerminationCriteria(const TerminationCriteria& criteria) { if (criteria.optimality_norm() != OPTIMALITY_NORM_L_INF && criteria.optimality_norm() != OPTIMALITY_NORM_L2) { return InvalidArgumentError("invalid value for optimality_norm"); } - if (std::isnan(criteria.eps_optimal_absolute())) { - return InvalidArgumentError("eps_optimal_absolute is NAN"); - } - if (criteria.eps_optimal_absolute() < 0) { - return InvalidArgumentError("eps_optimal_absolute must be non-negative"); - } - if (std::isnan(criteria.eps_optimal_relative())) { - return InvalidArgumentError("eps_optimal_relative is NAN"); - } - if (criteria.eps_optimal_relative() < 0) { - return InvalidArgumentError("eps_optimal_relative must be non-negative"); - } - if (std::isnan(criteria.eps_primal_infeasible())) { - return InvalidArgumentError("eps_primal_infeasible is NAN"); - } - if (criteria.eps_primal_infeasible() < 0) { - return InvalidArgumentError("eps_primal_infeasible must be non-negative"); - } - if (std::isnan(criteria.eps_dual_infeasible())) { - return InvalidArgumentError("eps_dual_infeasible is NAN"); - } - if (criteria.eps_dual_infeasible() < 0) { - return InvalidArgumentError("eps_dual_infeasible must be non-negative"); - } - if (std::isnan(criteria.time_sec_limit())) { - return InvalidArgumentError("time_sec_limit is NAN"); - } - if (criteria.time_sec_limit() < 0) { - return InvalidArgumentError("time_sec_limit must be non-negative"); - } + if (criteria.has_detailed_optimality_criteria() || + criteria.has_simple_optimality_criteria()) { + if (criteria.has_eps_optimal_absolute()) { + return InvalidArgumentError( + "eps_optimal_absolute should not be set if " + "detailed_optimality_criteria or simple_optimality_criteria is used"); + } + if (criteria.has_eps_optimal_relative()) { + return InvalidArgumentError( + "eps_optimal_relative should not be set if " + "detailed_optimality_criteria or simple_optimality_criteria is used"); + } + } + if (criteria.has_detailed_optimality_criteria()) { + RETURN_IF_ERROR(CheckNonNegative( + criteria.detailed_optimality_criteria() + .eps_optimal_primal_residual_absolute(), + "detailed_optimality_criteria.eps_optimal_primal_residual_absolute")); + RETURN_IF_ERROR(CheckNonNegative( + criteria.detailed_optimality_criteria() + .eps_optimal_primal_residual_relative(), + "detailed_optimality_criteria.eps_optimal_primal_residual_relative")); + RETURN_IF_ERROR(CheckNonNegative( + criteria.detailed_optimality_criteria() + .eps_optimal_dual_residual_absolute(), + "detailed_optimality_criteria.eps_optimal_dual_residual_absolute")); + RETURN_IF_ERROR(CheckNonNegative( + criteria.detailed_optimality_criteria() + .eps_optimal_dual_residual_relative(), + "detailed_optimality_criteria.eps_optimal_dual_residual_relative")); + RETURN_IF_ERROR(CheckNonNegative( + criteria.detailed_optimality_criteria() + .eps_optimal_objective_gap_absolute(), + "detailed_optimality_criteria.eps_optimal_objective_gap_absolute")); + RETURN_IF_ERROR(CheckNonNegative( + criteria.detailed_optimality_criteria() + .eps_optimal_objective_gap_relative(), + "detailed_optimality_criteria.eps_optimal_objective_gap_relative")); + } else if (criteria.has_simple_optimality_criteria()) { + RETURN_IF_ERROR(CheckNonNegative( + criteria.simple_optimality_criteria().eps_optimal_absolute(), + "simple_optimality_criteria.eps_optimal_absolute")); + RETURN_IF_ERROR(CheckNonNegative( + criteria.simple_optimality_criteria().eps_optimal_relative(), + "simple_optimality_criteria.eps_optimal_relative")); + } else { + RETURN_IF_ERROR(CheckNonNegative(criteria.eps_optimal_absolute(), + "eps_optimal_absolute")); + RETURN_IF_ERROR(CheckNonNegative(criteria.eps_optimal_relative(), + "eps_optimal_relative")); + } + RETURN_IF_ERROR(CheckNonNegative(criteria.eps_primal_infeasible(), + "eps_primal_infeasible")); + RETURN_IF_ERROR( + CheckNonNegative(criteria.eps_dual_infeasible(), "eps_dual_infeasible")); + RETURN_IF_ERROR( + CheckNonNegative(criteria.eps_dual_infeasible(), "eps_dual_infeasible")); + RETURN_IF_ERROR( + CheckNonNegative(criteria.time_sec_limit(), "time_sec_limit")); if (criteria.iteration_limit() < 0) { return InvalidArgumentError("iteration_limit must be non-negative"); } - if (std::isnan(criteria.kkt_matrix_pass_limit())) { - return InvalidArgumentError("kkt_matrix_pass_limit is NAN"); - } - if (criteria.kkt_matrix_pass_limit() < 0) { - return InvalidArgumentError("kkt_matrix_pass_limit must be non-negative"); - } + RETURN_IF_ERROR(CheckNonNegative(criteria.kkt_matrix_pass_limit(), + "kkt_matrix_pass_limit")); return OkStatus(); } diff --git a/ortools/pdlp/solvers_proto_validation_test.cc b/ortools/pdlp/solvers_proto_validation_test.cc index 96352ce3a43..039d358cb8f 100644 --- a/ortools/pdlp/solvers_proto_validation_test.cc +++ b/ortools/pdlp/solvers_proto_validation_test.cc @@ -17,13 +17,17 @@ #include #include "absl/status/status.h" +#include "absl/strings/str_cat.h" +//#include "absl/strings/string_view.h" #include "gmock/gmock.h" #include "gtest/gtest.h" +#include "ortools/base/protobuf_util.h" #include "ortools/pdlp/solvers.pb.h" namespace operations_research::pdlp { namespace { +using ::google::protobuf::util::ParseTextOrDie; using ::testing::HasSubstr; TEST(ValidateTerminationCriteria, DefaultIsValid) { @@ -40,36 +44,128 @@ TEST(ValidateTerminationCriteria, BadOptimalityNorm) { EXPECT_THAT(status.message(), HasSubstr("optimality_norm")); } -TEST(ValidateTerminationCriteria, BadEpsOptimalAbsolute) { - TerminationCriteria criteria_negative; - criteria_negative.set_eps_optimal_absolute(-1.0); - const absl::Status status_negative = - ValidateTerminationCriteria(criteria_negative); - EXPECT_EQ(status_negative.code(), absl::StatusCode::kInvalidArgument); - EXPECT_THAT(status_negative.message(), HasSubstr("eps_optimal_absolute")); +void TestTerminationCriteriaValidation( + std::string_view termination_criteria_string, + std::string_view error_substring) { + TerminationCriteria termination_criteria = + ParseTextOrDie(std::string(termination_criteria_string)); + const absl::Status status = ValidateTerminationCriteria(termination_criteria); + EXPECT_EQ(status.code(), absl::StatusCode::kInvalidArgument) + << "With termination criteria \"" << termination_criteria_string << "\""; + EXPECT_THAT(status.message(), HasSubstr(error_substring)) + << "With termination criteria \"" << termination_criteria_string << "\""; +} + +// Tests that the given DetailedOptimalityCriteria field can't be negative or +// NAN. +void TestDetailedOptimalityCriteriaFieldValidation( + std::string_view field_name) { + const std::string full_field_name = + absl::StrCat("detailed_optimality_criteria.", field_name); + TestTerminationCriteriaValidation( + absl::StrCat("detailed_optimality_criteria { ", field_name, ": -1.0 }"), + full_field_name); + TestTerminationCriteriaValidation( + absl::StrCat("detailed_optimality_criteria { ", field_name, ": nan }"), + full_field_name); +} - TerminationCriteria criteria_nan; - criteria_nan.set_eps_optimal_absolute( - std::numeric_limits::quiet_NaN()); - const absl::Status status_nan = ValidateTerminationCriteria(criteria_nan); - EXPECT_EQ(status_nan.code(), absl::StatusCode::kInvalidArgument); - EXPECT_THAT(status_nan.message(), HasSubstr("eps_optimal_absolute")); +TEST(ValidateTerminationCriteria, BadEpsOptimalAbsolute) { + TestTerminationCriteriaValidation("eps_optimal_absolute: -1.0", + "eps_optimal_absolute"); + TestTerminationCriteriaValidation("eps_optimal_absolute: nan", + "eps_optimal_absolute"); } TEST(ValidateTerminationCriteria, BadEpsOptimalRelative) { - TerminationCriteria criteria_negative; - criteria_negative.set_eps_optimal_relative(-1.0); - const absl::Status status_negative = - ValidateTerminationCriteria(criteria_negative); - EXPECT_EQ(status_negative.code(), absl::StatusCode::kInvalidArgument); - EXPECT_THAT(status_negative.message(), HasSubstr("eps_optimal_relative")); + TestTerminationCriteriaValidation("eps_optimal_relative: -1.0", + "eps_optimal_relative"); + TestTerminationCriteriaValidation("eps_optimal_relative: nan", + "eps_optimal_relative"); +} - TerminationCriteria criteria_nan; - criteria_nan.set_eps_optimal_relative( - std::numeric_limits::quiet_NaN()); - const absl::Status status_nan = ValidateTerminationCriteria(criteria_nan); - EXPECT_EQ(status_nan.code(), absl::StatusCode::kInvalidArgument); - EXPECT_THAT(status_nan.message(), HasSubstr("eps_optimal_relative")); +TEST(ValidateTerminationCriteria, BadSimpleEpsOptimalAbsolute) { + TestTerminationCriteriaValidation( + "simple_optimality_criteria { eps_optimal_absolute: -1.0}", + "simple_optimality_criteria.eps_optimal_absolute"); + TestTerminationCriteriaValidation( + "simple_optimality_criteria { eps_optimal_absolute: nan}", + "simple_optimality_criteria.eps_optimal_absolute"); +} + +TEST(ValidateTerminationCriteria, BadSimpleEpsOptimalRelative) { + TestTerminationCriteriaValidation( + "simple_optimality_criteria { eps_optimal_relative: -1.0}", + "simple_optimality_criteria.eps_optimal_relative"); + TestTerminationCriteriaValidation( + "simple_optimality_criteria { eps_optimal_relative: nan}", + "simple_optimality_criteria.eps_optimal_relative"); +} + +TEST(ValidateTerminationCriteria, BadDetailedEpsOptimalPrimalResidualAbsolute) { + TestDetailedOptimalityCriteriaFieldValidation( + "eps_optimal_primal_residual_absolute"); +} + +TEST(ValidateTerminationCriteria, BadDetailedEpsOptimalPrimalResidualRelative) { + TestDetailedOptimalityCriteriaFieldValidation( + "eps_optimal_primal_residual_relative"); +} + +TEST(ValidateTerminationCriteria, BadDetailedEpsOptimalDualResidualAbsolute) { + TestDetailedOptimalityCriteriaFieldValidation( + "eps_optimal_dual_residual_absolute"); +} + +TEST(ValidateTerminationCriteria, BadDetailedEpsOptimalDualResidualRelative) { + TestDetailedOptimalityCriteriaFieldValidation( + "eps_optimal_dual_residual_relative"); +} + +TEST(ValidateTerminationCriteria, BadDetailedEpsOptimalDualityGapAbsolute) { + TestDetailedOptimalityCriteriaFieldValidation( + "eps_optimal_objective_gap_absolute"); +} + +TEST(ValidateTerminationCriteria, BadDetailedEpsOptimalDualityGapRelative) { + TestDetailedOptimalityCriteriaFieldValidation( + "eps_optimal_objective_gap_relative"); +} + +TEST(ValidateTerminationCriteria, AbsoluteAndSimpleOptimalityCriteria) { + TerminationCriteria termination_criteria = + ParseTextOrDie( + "eps_optimal_absolute: 1.0 simple_optimality_criteria { }"); + const absl::Status status = ValidateTerminationCriteria(termination_criteria); + EXPECT_EQ(status.code(), absl::StatusCode::kInvalidArgument); + EXPECT_THAT(status.message(), HasSubstr("simple_optimality_criteria")); +} + +TEST(ValidateTerminationCriteria, RelativeAndSimpleOptimalityCriteria) { + TerminationCriteria termination_criteria = + ParseTextOrDie( + "eps_optimal_relative: 1.0 simple_optimality_criteria { }"); + const absl::Status status = ValidateTerminationCriteria(termination_criteria); + EXPECT_EQ(status.code(), absl::StatusCode::kInvalidArgument); + EXPECT_THAT(status.message(), HasSubstr("simple_optimality_criteria")); +} + +TEST(ValidateTerminationCriteria, AbsoluteAndDetailedOptimalityCriteria) { + TerminationCriteria termination_criteria = + ParseTextOrDie( + "eps_optimal_absolute: 1.0 detailed_optimality_criteria { }"); + const absl::Status status = ValidateTerminationCriteria(termination_criteria); + EXPECT_EQ(status.code(), absl::StatusCode::kInvalidArgument); + EXPECT_THAT(status.message(), HasSubstr("detailed_optimality_criteria")); +} + +TEST(ValidateTerminationCriteria, RelativeAndDetailedOptimalityCriteria) { + TerminationCriteria termination_criteria = + ParseTextOrDie( + "eps_optimal_relative: 1.0 detailed_optimality_criteria { }"); + const absl::Status status = ValidateTerminationCriteria(termination_criteria); + EXPECT_EQ(status.code(), absl::StatusCode::kInvalidArgument); + EXPECT_THAT(status.message(), HasSubstr("detailed_optimality_criteria")); } TEST(ValidateTerminationCriteria, BadEpsPriamlInfeasible) { diff --git a/ortools/pdlp/termination.cc b/ortools/pdlp/termination.cc index 010d76a0c41..4ff7a6160fd 100644 --- a/ortools/pdlp/termination.cc +++ b/ortools/pdlp/termination.cc @@ -13,6 +13,7 @@ #include "ortools/pdlp/termination.h" +#include #include #include #include @@ -25,10 +26,11 @@ namespace operations_research::pdlp { namespace { -bool OptimalityCriteriaMet(const OptimalityNorm optimality_norm, - const double abs_tol, const double rel_tol, - const ConvergenceInformation& stats, - const QuadraticProgramBoundNorms& bound_norms) { +bool OptimalityCriteriaMet( + const OptimalityNorm optimality_norm, + const TerminationCriteria::DetailedOptimalityCriteria& optimality_criteria, + const ConvergenceInformation& stats, + const QuadraticProgramBoundNorms& bound_norms) { const double abs_obj = std::abs(stats.primal_objective()) + std::abs(stats.dual_objective()); const double gap = @@ -57,9 +59,18 @@ bool OptimalityCriteriaMet(const OptimalityNorm optimality_norm, << OptimalityNorm_Name(optimality_norm); } - return dual_err <= abs_tol + rel_tol * dual_err_baseline && - primal_err <= abs_tol + rel_tol * primal_err_baseline && - std::isfinite(abs_obj) && gap <= abs_tol + rel_tol * abs_obj; + return dual_err <= + optimality_criteria.eps_optimal_dual_residual_absolute() + + optimality_criteria.eps_optimal_dual_residual_relative() * + dual_err_baseline && + primal_err <= + optimality_criteria.eps_optimal_primal_residual_absolute() + + optimality_criteria.eps_optimal_primal_residual_relative() * + primal_err_baseline && + std::isfinite(abs_obj) && + gap <= optimality_criteria.eps_optimal_objective_gap_absolute() + + optimality_criteria.eps_optimal_objective_gap_relative() * + abs_obj; } // Checks if the criteria for primal infeasibility are approximately @@ -86,6 +97,40 @@ bool DualInfeasibilityCriteriaMet(double eps_dual_infeasible, } } // namespace +TerminationCriteria::DetailedOptimalityCriteria EffectiveOptimalityCriteria( + const TerminationCriteria& termination_criteria) { + if (termination_criteria.has_detailed_optimality_criteria()) { + return termination_criteria.detailed_optimality_criteria(); + } + TerminationCriteria::SimpleOptimalityCriteria simple_criteria; + if (termination_criteria.has_simple_optimality_criteria()) { + simple_criteria = termination_criteria.simple_optimality_criteria(); + } else { + simple_criteria.set_eps_optimal_absolute( + termination_criteria.eps_optimal_absolute()); + simple_criteria.set_eps_optimal_relative( + termination_criteria.eps_optimal_relative()); + } + return EffectiveOptimalityCriteria(simple_criteria); +} + +TerminationCriteria::DetailedOptimalityCriteria EffectiveOptimalityCriteria( + const TerminationCriteria::SimpleOptimalityCriteria& simple_criteria) { + TerminationCriteria::DetailedOptimalityCriteria result; + result.set_eps_optimal_primal_residual_absolute( + simple_criteria.eps_optimal_absolute()); + result.set_eps_optimal_primal_residual_relative( + simple_criteria.eps_optimal_relative()); + result.set_eps_optimal_dual_residual_absolute( + simple_criteria.eps_optimal_absolute()); + result.set_eps_optimal_dual_residual_relative( + simple_criteria.eps_optimal_relative()); + result.set_eps_optimal_objective_gap_absolute( + simple_criteria.eps_optimal_absolute()); + result.set_eps_optimal_objective_gap_relative( + simple_criteria.eps_optimal_relative()); + return result; +} std::optional CheckSimpleTerminationCriteria( const TerminationCriteria& criteria, const IterationStats& stats, @@ -117,10 +162,11 @@ std::optional CheckTerminationCriteria( const QuadraticProgramBoundNorms& bound_norms, const std::atomic* interrupt_solve, const bool force_numerical_termination) { + TerminationCriteria::DetailedOptimalityCriteria optimality_criteria = + EffectiveOptimalityCriteria(criteria); for (const auto& convergence_stats : stats.convergence_information()) { - if (OptimalityCriteriaMet( - criteria.optimality_norm(), criteria.eps_optimal_absolute(), - criteria.eps_optimal_relative(), convergence_stats, bound_norms)) { + if (OptimalityCriteriaMet(criteria.optimality_norm(), optimality_criteria, + convergence_stats, bound_norms)) { return TerminationReasonAndPointType{ .reason = TERMINATION_REASON_OPTIMAL, .type = convergence_stats.candidate_type()}; @@ -149,7 +195,7 @@ std::optional CheckTerminationCriteria( return TerminationReasonAndPointType{ .reason = TERMINATION_REASON_NUMERICAL_ERROR, .type = POINT_TYPE_NONE}; } - return absl::nullopt; + return std::nullopt; } QuadraticProgramBoundNorms BoundNormsFromProblemStats( @@ -162,29 +208,39 @@ QuadraticProgramBoundNorms BoundNormsFromProblemStats( } RelativeConvergenceInformation ComputeRelativeResiduals( - const double eps_optimal_absolute, const double eps_optimal_relative, + const TerminationCriteria::DetailedOptimalityCriteria& optimality_criteria, const QuadraticProgramBoundNorms& norms, const ConvergenceInformation& stats) { - const double eps_ratio = eps_optimal_relative == 0.0 - ? std::numeric_limits::infinity() - : eps_optimal_absolute / eps_optimal_relative; + auto eps_ratio = [](const double absolute, const double relative) { + return relative == 0.0 ? std::numeric_limits::infinity() + : absolute / relative; + }; + const double eps_ratio_primal = + eps_ratio(optimality_criteria.eps_optimal_primal_residual_absolute(), + optimality_criteria.eps_optimal_primal_residual_relative()); + const double eps_ratio_dual = + eps_ratio(optimality_criteria.eps_optimal_dual_residual_absolute(), + optimality_criteria.eps_optimal_dual_residual_relative()); + const double eps_ratio_gap = + eps_ratio(optimality_criteria.eps_optimal_objective_gap_absolute(), + optimality_criteria.eps_optimal_objective_gap_relative()); RelativeConvergenceInformation info; info.relative_l_inf_primal_residual = stats.l_inf_primal_residual() / - (eps_ratio + norms.l_inf_norm_constraint_bounds); + (eps_ratio_primal + norms.l_inf_norm_constraint_bounds); info.relative_l2_primal_residual = stats.l2_primal_residual() / - (eps_ratio + norms.l2_norm_constraint_bounds); + (eps_ratio_primal + norms.l2_norm_constraint_bounds); info.relative_l_inf_dual_residual = stats.l_inf_dual_residual() / - (eps_ratio + norms.l_inf_norm_primal_linear_objective); + (eps_ratio_dual + norms.l_inf_norm_primal_linear_objective); info.relative_l2_dual_residual = stats.l2_dual_residual() / - (eps_ratio + norms.l2_norm_primal_linear_objective); + (eps_ratio_dual + norms.l2_norm_primal_linear_objective); const double abs_obj = std::abs(stats.primal_objective()) + std::abs(stats.dual_objective()); const double gap = stats.primal_objective() - stats.dual_objective(); - info.relative_optimality_gap = gap / (eps_ratio + abs_obj); + info.relative_optimality_gap = gap / (eps_ratio_gap + abs_obj); return info; } diff --git a/ortools/pdlp/termination.h b/ortools/pdlp/termination.h index ba302307fb2..af757c57383 100644 --- a/ortools/pdlp/termination.h +++ b/ortools/pdlp/termination.h @@ -14,6 +14,7 @@ #ifndef PDLP_TERMINATION_H_ #define PDLP_TERMINATION_H_ +#include #include #include "ortools/pdlp/solve_log.pb.h" @@ -35,6 +36,14 @@ struct QuadraticProgramBoundNorms { double l_inf_norm_constraint_bounds; }; +// Computes the effective optimality criteria for a TerminationCriteria. +TerminationCriteria::DetailedOptimalityCriteria EffectiveOptimalityCriteria( + const TerminationCriteria& termination_criteria); + +// Like previous overload but takes a SimpleOptimalityCriteria. Useful in +// unit tests where no TerminationCriteria is naturally available. +TerminationCriteria::DetailedOptimalityCriteria EffectiveOptimalityCriteria( + const TerminationCriteria::SimpleOptimalityCriteria& simple_criteria); // Checks if any of the simple termination criteria are satisfied by `stats`, // and returns a termination reason if so, and nullopt otherwise. The "simple" // termination criteria are `time_sec_limit`, `iteration_limit`, @@ -98,7 +107,7 @@ struct RelativeConvergenceInformation { }; RelativeConvergenceInformation ComputeRelativeResiduals( - double eps_optimal_absolute, double eps_optimal_relative, + const TerminationCriteria::DetailedOptimalityCriteria& optimality_criteria, const QuadraticProgramBoundNorms& norms, const ConvergenceInformation& stats); diff --git a/ortools/pdlp/termination_test.cc b/ortools/pdlp/termination_test.cc index ffda5327bcc..b0115fbbb58 100644 --- a/ortools/pdlp/termination_test.cc +++ b/ortools/pdlp/termination_test.cc @@ -13,6 +13,7 @@ #include "ortools/pdlp/termination.h" +#include #include #include @@ -26,6 +27,7 @@ namespace operations_research::pdlp { namespace { using ::google::protobuf::util::ParseTextOrDie; +using ::testing::EqualsProto; using ::testing::FieldsAre; using ::testing::Optional; @@ -52,8 +54,10 @@ class TerminationTest : public testing::TestWithParam { protected: void SetUp() override { test_criteria_ = ParseTextOrDie(R"pb( - eps_optimal_absolute: 1.0e-4 - eps_optimal_relative: 1.0e-4 + simple_optimality_criteria { + eps_optimal_absolute: 1.0e-4 + eps_optimal_relative: 1.0e-4 + } eps_primal_infeasible: 1.0e-6 eps_dual_infeasible: 1.0e-6 time_sec_limit: 1.0 @@ -64,6 +68,230 @@ class TerminationTest : public testing::TestWithParam { TerminationCriteria test_criteria_; }; +class DetailedRelativeTerminationTest + : public testing::TestWithParam { + protected: + void SetUp() override { + test_criteria_ = ParseTextOrDie(R"pb( + detailed_optimality_criteria { + eps_optimal_primal_residual_absolute: 0.0 + eps_optimal_primal_residual_relative: 1.0e-4 + eps_optimal_dual_residual_absolute: 0.0 + eps_optimal_dual_residual_relative: 1.0e-4 + eps_optimal_objective_gap_absolute: 0.0 + eps_optimal_objective_gap_relative: 1.0e-4 + } + )pb"); + test_criteria_.set_optimality_norm(GetParam()); + } + + TerminationCriteria test_criteria_; +}; + +class DetailedAbsoluteTerminationTest + : public testing::TestWithParam { + protected: + void SetUp() override { + test_criteria_ = ParseTextOrDie(R"pb( + detailed_optimality_criteria { + eps_optimal_primal_residual_absolute: 1.0e-4 + eps_optimal_primal_residual_relative: 0.0 + eps_optimal_dual_residual_absolute: 1.0e-4 + eps_optimal_dual_residual_relative: 0.0 + eps_optimal_objective_gap_absolute: 1.0e-4 + eps_optimal_objective_gap_relative: 0.0 + } + )pb"); + test_criteria_.set_optimality_norm(GetParam()); + } + + TerminationCriteria test_criteria_; +}; + +TEST(EffectiveOptimalityCriteriaTest, SimpleOptimalityCriteriaOverload) { + const auto criteria = + ParseTextOrDie( + R"pb(eps_optimal_absolute: 1.0e-4 eps_optimal_relative: 2.0e-4)pb"); + EXPECT_THAT(EffectiveOptimalityCriteria(criteria), EqualsProto(R"pb( + eps_optimal_primal_residual_absolute: 1.0e-4 + eps_optimal_primal_residual_relative: 2.0e-4 + eps_optimal_dual_residual_absolute: 1.0e-4 + eps_optimal_dual_residual_relative: 2.0e-4 + eps_optimal_objective_gap_absolute: 1.0e-4 + eps_optimal_objective_gap_relative: 2.0e-4 + )pb")); +} + +TEST(EffectiveOptimalityCriteriaTest, SimpleOptimalityCriteriaInput) { + const auto criteria = + ParseTextOrDie(R"pb(simple_optimality_criteria { + eps_optimal_absolute: 1.0e-4 + eps_optimal_relative: 2.0e-4 + })pb"); + EXPECT_THAT(EffectiveOptimalityCriteria(criteria), EqualsProto(R"pb( + eps_optimal_primal_residual_absolute: 1.0e-4 + eps_optimal_primal_residual_relative: 2.0e-4 + eps_optimal_dual_residual_absolute: 1.0e-4 + eps_optimal_dual_residual_relative: 2.0e-4 + eps_optimal_objective_gap_absolute: 1.0e-4 + eps_optimal_objective_gap_relative: 2.0e-4 + )pb")); +} + +TEST(EffectiveOptimalityCriteriaTest, DetailedOptimalityCriteriaInput) { + const auto criteria = ParseTextOrDie( + R"pb(detailed_optimality_criteria { + eps_optimal_primal_residual_absolute: 1.0e-4 + eps_optimal_primal_residual_relative: 2.0e-4 + eps_optimal_dual_residual_absolute: 3.0e-4 + eps_optimal_dual_residual_relative: 4.0e-4 + eps_optimal_objective_gap_absolute: 5.0e-4 + eps_optimal_objective_gap_relative: 6.0e-4 + })pb"); + EXPECT_THAT(EffectiveOptimalityCriteria(criteria), + EqualsProto(criteria.detailed_optimality_criteria())); +} + +TEST(EffectiveOptimalityCriteriaTest, DeprecatedInput) { + const auto criteria = ParseTextOrDie( + R"pb(eps_optimal_absolute: 1.0e-4 eps_optimal_relative: 2.0e-4)pb"); + EXPECT_THAT(EffectiveOptimalityCriteria(criteria), EqualsProto(R"pb( + eps_optimal_primal_residual_absolute: 1.0e-4 + eps_optimal_primal_residual_relative: 2.0e-4 + eps_optimal_dual_residual_absolute: 1.0e-4 + eps_optimal_dual_residual_relative: 2.0e-4 + eps_optimal_objective_gap_absolute: 1.0e-4 + eps_optimal_objective_gap_relative: 2.0e-4 + )pb")); +} + +TEST_P(DetailedRelativeTerminationTest, TerminationWithNearOptimal) { + const auto stats = ParseTextOrDie(R"pb( + convergence_information { + primal_objective: 1.00019 + dual_objective: 1.0 + l_inf_primal_residual: 11e-4 + l_inf_dual_residual: 5.4e-4 + l2_primal_residual: 14e-4 + l2_dual_residual: 6e-4 + candidate_type: POINT_TYPE_CURRENT_ITERATE + })pb"); + EXPECT_THAT( + CheckTerminationCriteria(test_criteria_, stats, TestLpBoundNorms()), + Optional( + FieldsAre(TERMINATION_REASON_OPTIMAL, POINT_TYPE_CURRENT_ITERATE))); +} + +TEST_P(DetailedRelativeTerminationTest, NoTerminationWithExcessiveGap) { + const auto stats = ParseTextOrDie(R"pb( + convergence_information { + primal_objective: 1.00021 + dual_objective: 1.0 + l_inf_primal_residual: 11e-4 + l_inf_dual_residual: 5.4e-4 + l2_primal_residual: 14e-4 + l2_dual_residual: 6e-4 + candidate_type: POINT_TYPE_CURRENT_ITERATE + })pb"); + EXPECT_EQ(CheckTerminationCriteria(test_criteria_, stats, TestLpBoundNorms()), + absl::nullopt); +} + +TEST_P(DetailedRelativeTerminationTest, + NoTerminationWithExcessivePrimalResidual) { + const auto stats = ParseTextOrDie(R"pb( + convergence_information { + primal_objective: 1.00019 + dual_objective: 1.0 + l_inf_primal_residual: 13e-4 + l_inf_dual_residual: 5.4e-4 + l2_primal_residual: 15e-4 + l2_dual_residual: 6e-4 + candidate_type: POINT_TYPE_CURRENT_ITERATE + })pb"); + EXPECT_EQ(CheckTerminationCriteria(test_criteria_, stats, TestLpBoundNorms()), + absl::nullopt); +} + +TEST_P(DetailedRelativeTerminationTest, + NoTerminationWithExcessiveDualResidual) { + const auto stats = ParseTextOrDie(R"pb( + convergence_information { + primal_objective: 1.00019 + dual_objective: 1.0 + l_inf_primal_residual: 11e-4 + l_inf_dual_residual: 5.6e-4 + l2_primal_residual: 14e-4 + l2_dual_residual: 7e-4 + candidate_type: POINT_TYPE_CURRENT_ITERATE + })pb"); + EXPECT_EQ(CheckTerminationCriteria(test_criteria_, stats, TestLpBoundNorms()), + absl::nullopt); +} + +TEST_P(DetailedAbsoluteTerminationTest, TerminationWithNearOptimal) { + const auto stats = ParseTextOrDie(R"pb( + convergence_information { + primal_objective: 1.00009 + dual_objective: 1.0 + l_inf_primal_residual: 9e-5 + l_inf_dual_residual: 9e-5 + l2_primal_residual: 9e-5 + l2_dual_residual: 9e-5 + candidate_type: POINT_TYPE_CURRENT_ITERATE + })pb"); + EXPECT_THAT( + CheckTerminationCriteria(test_criteria_, stats, TestLpBoundNorms()), + Optional( + FieldsAre(TERMINATION_REASON_OPTIMAL, POINT_TYPE_CURRENT_ITERATE))); +} + +TEST_P(DetailedAbsoluteTerminationTest, NoTerminationWithExcessiveGap) { + const auto stats = ParseTextOrDie(R"pb( + convergence_information { + primal_objective: 1.00011 + dual_objective: 1.0 + l_inf_primal_residual: 9e-5 + l_inf_dual_residual: 9e-5 + l2_primal_residual: 9e-5 + l2_dual_residual: 9e-5 + candidate_type: POINT_TYPE_CURRENT_ITERATE + })pb"); + EXPECT_EQ(CheckTerminationCriteria(test_criteria_, stats, TestLpBoundNorms()), + absl::nullopt); +} + +TEST_P(DetailedAbsoluteTerminationTest, + NoTerminationWithExcessivePrimalResidual) { + const auto stats = ParseTextOrDie(R"pb( + convergence_information { + primal_objective: 1.00009 + dual_objective: 1.0 + l_inf_primal_residual: 11e-5 + l_inf_dual_residual: 9e-5 + l2_primal_residual: 11e-5 + l2_dual_residual: 9e-5 + candidate_type: POINT_TYPE_CURRENT_ITERATE + })pb"); + EXPECT_EQ(CheckTerminationCriteria(test_criteria_, stats, TestLpBoundNorms()), + absl::nullopt); +} + +TEST_P(DetailedAbsoluteTerminationTest, + NoTerminationWithExcessiveDualResidual) { + const auto stats = ParseTextOrDie(R"pb( + convergence_information { + primal_objective: 1.00009 + dual_objective: 1.0 + l_inf_primal_residual: 9e-5 + l_inf_dual_residual: 11e-5 + l2_primal_residual: 9e-5 + l2_dual_residual: 11e-5 + candidate_type: POINT_TYPE_CURRENT_ITERATE + })pb"); + EXPECT_EQ(CheckTerminationCriteria(test_criteria_, stats, TestLpBoundNorms()), + absl::nullopt); +} TEST_P(TerminationTest, NoTerminationWithLargeGap) { IterationStats stats = ParseTextOrDie(R"pb( @@ -73,19 +301,19 @@ TEST_P(TerminationTest, NoTerminationWithLargeGap) { dual_objective: -50.0 })pb"); EXPECT_EQ(CheckTerminationCriteria(test_criteria_, stats, TestLpBoundNorms()), - absl::nullopt); + std::nullopt); } TEST_F(SimpleTerminationTest, NoTerminationWithEmptyIterationStats) { IterationStats stats; EXPECT_EQ(CheckSimpleTerminationCriteria(test_criteria_, stats), - absl::nullopt); + std::nullopt); } TEST_P(TerminationTest, NoTerminationWithEmptyIterationStats) { IterationStats stats; EXPECT_EQ(CheckTerminationCriteria(test_criteria_, stats, TestLpBoundNorms()), - absl::nullopt); + std::nullopt); } TEST_F(SimpleTerminationTest, TerminationWithInterruptSolve) { @@ -170,7 +398,7 @@ TEST_P(TerminationTest, NoTerminationWithInfeasibleDualRay) { })pb"); EXPECT_EQ(CheckTerminationCriteria(test_criteria_, stats_infeasible_ray, TestLpBoundNorms()), - absl::nullopt); + std::nullopt); } TEST_P(TerminationTest, NoTerminationWithNegativeDualRayObjective) { @@ -181,7 +409,7 @@ TEST_P(TerminationTest, NoTerminationWithNegativeDualRayObjective) { })pb"); EXPECT_EQ(CheckTerminationCriteria(test_criteria_, stats_wrong_sign, TestLpBoundNorms()), - absl::nullopt); + std::nullopt); } TEST_P(TerminationTest, NoTerminationWithZeroDualRayObjective) { @@ -192,7 +420,7 @@ TEST_P(TerminationTest, NoTerminationWithZeroDualRayObjective) { })pb"); EXPECT_EQ(CheckTerminationCriteria(test_criteria_, stats_objective_zero, TestLpBoundNorms()), - absl::nullopt); + std::nullopt); } TEST_P(TerminationTest, DualInfeasibleFromAverageIterate) { @@ -217,7 +445,7 @@ TEST_P(TerminationTest, NoTerminationWithInfeasiblePrimalRay) { })pb"); EXPECT_EQ(CheckTerminationCriteria(test_criteria_, stats_infeasible_ray, TestLpBoundNorms()), - absl::nullopt); + std::nullopt); } TEST_P(TerminationTest, NoTerminationWithPositivePrimalRayObjective) { @@ -228,7 +456,7 @@ TEST_P(TerminationTest, NoTerminationWithPositivePrimalRayObjective) { })pb"); EXPECT_EQ(CheckTerminationCriteria(test_criteria_, stats_wrong_sign, TestLpBoundNorms()), - absl::nullopt); + std::nullopt); } TEST_P(TerminationTest, NoTerminationWithZeroPrimalRayObjective) { @@ -239,10 +467,10 @@ TEST_P(TerminationTest, NoTerminationWithZeroPrimalRayObjective) { })pb"); EXPECT_EQ(CheckTerminationCriteria(test_criteria_, stats_objective_zero, TestLpBoundNorms()), - absl::nullopt); + std::nullopt); } -TEST_P(TerminationTest, Optimal) { +TEST_P(TerminationTest, TerminationWithOptimal) { const auto stats = ParseTextOrDie(R"pb( convergence_information { primal_objective: 1.0 @@ -260,6 +488,24 @@ TEST_P(TerminationTest, Optimal) { POINT_TYPE_CURRENT_ITERATE))); } +TEST_P(TerminationTest, TerminationWithNearOptimal) { + const auto stats = ParseTextOrDie(R"pb( + convergence_information { + primal_objective: 1.00019 + dual_objective: 1.0 + l_inf_primal_residual: 11e-4 + l_inf_dual_residual: 5.4e-4 + l2_primal_residual: 14e-4 + l2_dual_residual: 6e-4 + candidate_type: POINT_TYPE_CURRENT_ITERATE + })pb"); + + std::optional maybe_result = + CheckTerminationCriteria(test_criteria_, stats, TestLpBoundNorms()); + EXPECT_THAT(maybe_result, Optional(FieldsAre(TERMINATION_REASON_OPTIMAL, + POINT_TYPE_CURRENT_ITERATE))); +} + TEST_P(TerminationTest, OptimalEvenWithNumericalError) { const auto stats = ParseTextOrDie(R"pb( convergence_information { @@ -313,7 +559,7 @@ TEST_P(TerminationTest, NoTerminationWithBadGap) { })pb"); EXPECT_EQ(CheckTerminationCriteria(test_criteria_, stats_bad_gap, TestLpBoundNorms()), - absl::nullopt); + std::nullopt); } TEST_P(TerminationTest, NoTerminationWithInfiniteGap) { @@ -328,7 +574,7 @@ TEST_P(TerminationTest, NoTerminationWithInfiniteGap) { })pb"); EXPECT_EQ(CheckTerminationCriteria(test_criteria_, stats_infinite_gap, TestLpBoundNorms()), - absl::nullopt); + std::nullopt); } TEST_P(TerminationTest, NoTerminationWithBadPrimalResidual) { @@ -343,7 +589,7 @@ TEST_P(TerminationTest, NoTerminationWithBadPrimalResidual) { })pb"); EXPECT_EQ(CheckTerminationCriteria(test_criteria_, stats_bad_primal, TestLpBoundNorms()), - absl::nullopt); + std::nullopt); } TEST_P(TerminationTest, NoTerminationWithBadDualResidual) { @@ -358,7 +604,7 @@ TEST_P(TerminationTest, NoTerminationWithBadDualResidual) { })pb"); EXPECT_EQ(CheckTerminationCriteria(test_criteria_, stats_bad_dual, TestLpBoundNorms()), - absl::nullopt); + std::nullopt); } // Tests that optimality is checked with non-strict inequalities, as per the @@ -375,8 +621,10 @@ TEST_P(TerminationTest, ZeroToleranceZeroError) { candidate_type: POINT_TYPE_CURRENT_ITERATE })pb"); - test_criteria_.set_eps_optimal_absolute(0.0); - test_criteria_.set_eps_optimal_relative(0.0); + test_criteria_.mutable_simple_optimality_criteria()->set_eps_optimal_absolute( + 0.0); + test_criteria_.mutable_simple_optimality_criteria()->set_eps_optimal_relative( + 0.0); std::optional maybe_result = CheckTerminationCriteria(test_criteria_, stats, TestLpBoundNorms()); @@ -388,9 +636,18 @@ INSTANTIATE_TEST_SUITE_P(OptNorm, TerminationTest, testing::Values(OPTIMALITY_NORM_L2, OPTIMALITY_NORM_L_INF)); +INSTANTIATE_TEST_SUITE_P(DetailedRelativeOptNorm, + DetailedRelativeTerminationTest, + testing::Values(OPTIMALITY_NORM_L2, + OPTIMALITY_NORM_L_INF)); +INSTANTIATE_TEST_SUITE_P(DetailedAbsoluteOptNorm, + DetailedAbsoluteTerminationTest, + testing::Values(OPTIMALITY_NORM_L2, + OPTIMALITY_NORM_L_INF)); + TEST(TerminationTest, L2AndLInfDiffer) { auto test_criteria = ParseTextOrDie(R"pb( - eps_optimal_relative: 1.0)pb"); + simple_optimality_criteria { eps_optimal_relative: 1.0 })pb"); // For L2, optimality requires norm(primal_residual, 2) <= 14.49 // For LInf, optimality requires norm(primal_residual, Inf) <= 12.0 @@ -408,8 +665,8 @@ TEST(TerminationTest, L2AndLInfDiffer) { {13.0, TerminationReasonAndPointType{.reason = TERMINATION_REASON_OPTIMAL, .type = POINT_TYPE_CURRENT_ITERATE}, - absl::nullopt}, - {15.0, absl::nullopt, absl::nullopt}}; + std::nullopt}, + {15.0, std::nullopt, std::nullopt}}; for (const auto& config : test_configs) { IterationStats stats; @@ -469,9 +726,14 @@ TEST(ComputeRelativeResiduals, stats.set_l2_primal_residual(1.0); stats.set_l_inf_dual_residual(1.0); stats.set_l2_dual_residual(1.0); + TerminationCriteria termination_criteria; + termination_criteria.mutable_simple_optimality_criteria() + ->set_eps_optimal_absolute(0.0); + termination_criteria.mutable_simple_optimality_criteria() + ->set_eps_optimal_relative(1.0e-6); const RelativeConvergenceInformation relative_info = ComputeRelativeResiduals( - /*eps_optimal_absolute=*/0.0, /*eps_optimal_relative=*/1.0e-6, - TestLpBoundNorms(), stats); + EffectiveOptimalityCriteria(termination_criteria), TestLpBoundNorms(), + stats); EXPECT_EQ(relative_info.relative_l_inf_primal_residual, 1.0 / 12.0); EXPECT_EQ(relative_info.relative_l2_primal_residual, 1.0 / std::sqrt(210.0)); @@ -497,9 +759,11 @@ TEST(ComputeRelativeResiduals, stats.set_l2_primal_residual(1.0); stats.set_l_inf_dual_residual(1.0); stats.set_l2_dual_residual(1.0); + TerminationCriteria::SimpleOptimalityCriteria opt_criteria; + opt_criteria.set_eps_optimal_absolute(0.0); + opt_criteria.set_eps_optimal_relative(0.0); const RelativeConvergenceInformation relative_info = ComputeRelativeResiduals( - /*eps_optimal_absolute=*/0.0, /*eps_optimal_relative=*/0.0, - TestLpBoundNorms(), stats); + EffectiveOptimalityCriteria(opt_criteria), TestLpBoundNorms(), stats); EXPECT_EQ(relative_info.relative_l_inf_primal_residual, 0.0); EXPECT_EQ(relative_info.relative_l2_primal_residual, 0.0); @@ -520,9 +784,11 @@ TEST(ComputeRelativeResiduals, stats.set_l2_primal_residual(1.0); stats.set_l_inf_dual_residual(1.0); stats.set_l2_dual_residual(1.0); + TerminationCriteria::SimpleOptimalityCriteria opt_criteria; + opt_criteria.set_eps_optimal_absolute(1.0e-6); + opt_criteria.set_eps_optimal_relative(1.0e-6); const RelativeConvergenceInformation relative_info = ComputeRelativeResiduals( - /*eps_optimal_absolute=*/1.0e-6, /*eps_optimal_relative=*/1.0e-6, - TestLpBoundNorms(), stats); + EffectiveOptimalityCriteria(opt_criteria), TestLpBoundNorms(), stats); EXPECT_EQ(relative_info.relative_l_inf_primal_residual, 1.0 / (1.0 + 12.0)); EXPECT_EQ(relative_info.relative_l2_primal_residual, @@ -536,5 +802,40 @@ TEST(ComputeRelativeResiduals, EXPECT_EQ(relative_info.relative_optimality_gap, 5.0 / (1.0 + 15.0)); } +TEST(ComputeRelativeResiduals, + ComputesCorrectRelativeResidualsForDetailedTerminationCriteria) { + ConvergenceInformation stats; + // If the absolute error tolerance and relative error tolerance are equal (and + // nonzero), the relative residuals are the absolute residuals divided by 1.0 + // plus the corresponding norms. + stats.set_primal_objective(10.0); + stats.set_dual_objective(5.0); + stats.set_l_inf_primal_residual(1.0); + stats.set_l2_primal_residual(1.0); + stats.set_l_inf_dual_residual(1.0); + stats.set_l2_dual_residual(1.0); + TerminationCriteria::DetailedOptimalityCriteria opt_criteria; + opt_criteria.set_eps_optimal_primal_residual_absolute(2.0e-6); + opt_criteria.set_eps_optimal_primal_residual_relative(2.0e-4); + opt_criteria.set_eps_optimal_dual_residual_absolute(1.0e-3); + opt_criteria.set_eps_optimal_dual_residual_relative(1.0e-4); + opt_criteria.set_eps_optimal_objective_gap_absolute(3.0e-8); + opt_criteria.set_eps_optimal_objective_gap_relative(3.0e-7); + const RelativeConvergenceInformation relative_info = + ComputeRelativeResiduals(opt_criteria, TestLpBoundNorms(), stats); + + EXPECT_EQ(relative_info.relative_l_inf_primal_residual, 1.0 / (0.01 + 12.0)); + EXPECT_EQ(relative_info.relative_l2_primal_residual, + 1.0 / (0.01 + std::sqrt(210.0))); + + EXPECT_EQ(relative_info.relative_l_inf_dual_residual, 1.0 / (10.0 + 5.5)); + EXPECT_EQ(relative_info.relative_l2_dual_residual, + 1.0 / (10.0 + sqrt(36.25))); + + // The relative optimality gap should just be the objective difference divided + // by 0.1 + the sum of absolute values. + EXPECT_EQ(relative_info.relative_optimality_gap, 5.0 / (0.1 + 15.0)); +} + } // namespace } // namespace operations_research::pdlp diff --git a/ortools/pdlp/test_util.cc b/ortools/pdlp/test_util.cc index 36aebdcd2f9..c78c49dce4d 100644 --- a/ortools/pdlp/test_util.cc +++ b/ortools/pdlp/test_util.cc @@ -21,6 +21,7 @@ #include "Eigen/Core" #include "Eigen/SparseCore" #include "gmock/gmock.h" +#include "gtest/gtest.h" #include "ortools/pdlp/quadratic_program.h" namespace operations_research::pdlp { diff --git a/ortools/pdlp/test_util.h b/ortools/pdlp/test_util.h index 2a723d8df39..583ec028688 100644 --- a/ortools/pdlp/test_util.h +++ b/ortools/pdlp/test_util.h @@ -17,14 +17,14 @@ #include #include #include -#include +#include #include "Eigen/Core" #include "Eigen/SparseCore" #include "absl/types/span.h" #include "gmock/gmock.h" #include "gtest/gtest.h" -#include "ortools/base/logging.h" +#include "ortools/base/check.h" #include "ortools/pdlp/quadratic_program.h" namespace operations_research::pdlp { diff --git a/ortools/pdlp/test_util_test.cc b/ortools/pdlp/test_util_test.cc index 46335a78b40..7890b59989c 100644 --- a/ortools/pdlp/test_util_test.cc +++ b/ortools/pdlp/test_util_test.cc @@ -23,7 +23,7 @@ #include "absl/types/span.h" #include "gmock/gmock.h" #include "gtest/gtest.h" -#include "ortools/base/logging.h" +#include "ortools/base/check.h" namespace operations_research::pdlp { namespace { diff --git a/ortools/pdlp/trust_region.cc b/ortools/pdlp/trust_region.cc index ed1d96eabb5..8306f78eeb9 100644 --- a/ortools/pdlp/trust_region.cc +++ b/ortools/pdlp/trust_region.cc @@ -22,6 +22,7 @@ #include #include "Eigen/Core" +#include "ortools/base/check.h" #include "ortools/base/logging.h" #include "ortools/base/mathutil.h" #include "ortools/pdlp/quadratic_program.h" @@ -179,7 +180,7 @@ double MedianOfShardMedians( const std::vector>& indexed_components_by_shard, const Sharder& sharder) { std::vector> shard_medians(sharder.NumShards(), - absl::nullopt); + std::nullopt); sharder.ParallelForEachShard([&](const Sharder::Shard& shard) { const auto& indexed_shard_components = indexed_components_by_shard[shard.Index()]; diff --git a/ortools/pdlp/trust_region.h b/ortools/pdlp/trust_region.h index 9593a8c12c3..06123b208cc 100644 --- a/ortools/pdlp/trust_region.h +++ b/ortools/pdlp/trust_region.h @@ -22,7 +22,7 @@ #include "Eigen/Core" #include "absl/algorithm/container.h" -#include "ortools/base/logging.h" +#include "ortools/base/check.h" #include "ortools/base/mathutil.h" #include "ortools/pdlp/sharded_quadratic_program.h" #include "ortools/pdlp/sharder.h" diff --git a/ortools/pdlp/trust_region_test.cc b/ortools/pdlp/trust_region_test.cc index b78663c6ae6..2d0a83aca9a 100644 --- a/ortools/pdlp/trust_region_test.cc +++ b/ortools/pdlp/trust_region_test.cc @@ -17,10 +17,13 @@ #include #include #include +#include #include #include "Eigen/Core" #include "Eigen/SparseCore" +#include "absl/strings/str_cat.h" +#include "absl/strings/string_view.h" #include "gmock/gmock.h" #include "gtest/gtest.h" #include "ortools/pdlp/quadratic_program.h" From a1e4b6388f4130df0b4fd88d9fc245cbdbdfd691 Mon Sep 17 00:00:00 2001 From: Mizux Seiha Date: Tue, 20 Sep 2022 13:06:55 +0200 Subject: [PATCH 2/2] Fix pdlp export --- ortools/pdlp/termination_test.cc | 50 +++++++++++++++++++++++++++----- 1 file changed, 42 insertions(+), 8 deletions(-) diff --git a/ortools/pdlp/termination_test.cc b/ortools/pdlp/termination_test.cc index b0115fbbb58..c71a34f68e6 100644 --- a/ortools/pdlp/termination_test.cc +++ b/ortools/pdlp/termination_test.cc @@ -24,10 +24,41 @@ #include "ortools/pdlp/solvers.pb.h" namespace operations_research::pdlp { + +bool operator==( + const TerminationCriteria::DetailedOptimalityCriteria& lhs, + const TerminationCriteria::DetailedOptimalityCriteria& rhs) { + if (lhs.eps_optimal_primal_residual_absolute() != + rhs.eps_optimal_primal_residual_absolute()) { + return false; + } + if (lhs.eps_optimal_primal_residual_relative() != + rhs.eps_optimal_primal_residual_relative()) { + return false; + } + if (lhs.eps_optimal_dual_residual_absolute() != + rhs.eps_optimal_dual_residual_absolute()) { + return false; + } + if (lhs.eps_optimal_dual_residual_relative() != + rhs.eps_optimal_dual_residual_relative()) { + return false; + } + if (lhs.eps_optimal_objective_gap_absolute() != + rhs.eps_optimal_objective_gap_absolute()) { + return false; + } + if (lhs.eps_optimal_objective_gap_relative() != + rhs.eps_optimal_objective_gap_relative()) { + return false; + } + return true; +} + namespace { using ::google::protobuf::util::ParseTextOrDie; -using ::testing::EqualsProto; +using ::testing::Eq; using ::testing::FieldsAre; using ::testing::Optional; @@ -112,14 +143,15 @@ TEST(EffectiveOptimalityCriteriaTest, SimpleOptimalityCriteriaOverload) { const auto criteria = ParseTextOrDie( R"pb(eps_optimal_absolute: 1.0e-4 eps_optimal_relative: 2.0e-4)pb"); - EXPECT_THAT(EffectiveOptimalityCriteria(criteria), EqualsProto(R"pb( + EXPECT_THAT(EffectiveOptimalityCriteria(criteria), + Eq(ParseTextOrDie(R"pb( eps_optimal_primal_residual_absolute: 1.0e-4 eps_optimal_primal_residual_relative: 2.0e-4 eps_optimal_dual_residual_absolute: 1.0e-4 eps_optimal_dual_residual_relative: 2.0e-4 eps_optimal_objective_gap_absolute: 1.0e-4 eps_optimal_objective_gap_relative: 2.0e-4 - )pb")); + )pb"))); } TEST(EffectiveOptimalityCriteriaTest, SimpleOptimalityCriteriaInput) { @@ -128,14 +160,15 @@ TEST(EffectiveOptimalityCriteriaTest, SimpleOptimalityCriteriaInput) { eps_optimal_absolute: 1.0e-4 eps_optimal_relative: 2.0e-4 })pb"); - EXPECT_THAT(EffectiveOptimalityCriteria(criteria), EqualsProto(R"pb( + EXPECT_THAT(EffectiveOptimalityCriteria(criteria), + Eq(ParseTextOrDie(R"pb( eps_optimal_primal_residual_absolute: 1.0e-4 eps_optimal_primal_residual_relative: 2.0e-4 eps_optimal_dual_residual_absolute: 1.0e-4 eps_optimal_dual_residual_relative: 2.0e-4 eps_optimal_objective_gap_absolute: 1.0e-4 eps_optimal_objective_gap_relative: 2.0e-4 - )pb")); + )pb"))); } TEST(EffectiveOptimalityCriteriaTest, DetailedOptimalityCriteriaInput) { @@ -149,20 +182,21 @@ TEST(EffectiveOptimalityCriteriaTest, DetailedOptimalityCriteriaInput) { eps_optimal_objective_gap_relative: 6.0e-4 })pb"); EXPECT_THAT(EffectiveOptimalityCriteria(criteria), - EqualsProto(criteria.detailed_optimality_criteria())); + Eq(criteria.detailed_optimality_criteria())); } TEST(EffectiveOptimalityCriteriaTest, DeprecatedInput) { const auto criteria = ParseTextOrDie( R"pb(eps_optimal_absolute: 1.0e-4 eps_optimal_relative: 2.0e-4)pb"); - EXPECT_THAT(EffectiveOptimalityCriteria(criteria), EqualsProto(R"pb( + EXPECT_THAT(EffectiveOptimalityCriteria(criteria), + Eq(ParseTextOrDie(R"pb( eps_optimal_primal_residual_absolute: 1.0e-4 eps_optimal_primal_residual_relative: 2.0e-4 eps_optimal_dual_residual_absolute: 1.0e-4 eps_optimal_dual_residual_relative: 2.0e-4 eps_optimal_objective_gap_absolute: 1.0e-4 eps_optimal_objective_gap_relative: 2.0e-4 - )pb")); + )pb"))); } TEST_P(DetailedRelativeTerminationTest, TerminationWithNearOptimal) {