diff --git a/cpp/src/dual_simplex/presolve.cpp b/cpp/src/dual_simplex/presolve.cpp index e87c88b40..d94fd0f6c 100644 --- a/cpp/src/dual_simplex/presolve.cpp +++ b/cpp/src/dual_simplex/presolve.cpp @@ -842,6 +842,30 @@ void uncrush_primal_solution(const user_problem_t& user_problem, std::copy(solution.begin(), solution.begin() + user_problem.num_cols, user_solution.data()); } +template +void uncrush_dual_solution(const user_problem_t& user_problem, + const lp_problem_t& problem, + const std::vector& y, + const std::vector& z, + std::vector& user_y, + std::vector& user_z) +{ + // Reduced costs are uncrushed just like the primal solution + uncrush_primal_solution(user_problem, problem, z, user_z); + + // Adjust the sign of the dual variables y + // We should have A^T y + z = c + // In convert_user_problem, we converted >= to <=, so we need to adjust the sign of the dual + // variables + for (i_t i = 0; i < user_problem.num_rows; i++) { + if (user_problem.row_sense[i] == 'G') { + user_y[i] = -y[i]; + } else { + user_y[i] = y[i]; + } + } +} + template void uncrush_solution(const presolve_info_t& presolve_info, const std::vector& crushed_x, @@ -903,6 +927,13 @@ template void uncrush_primal_solution(const user_problem_t& solution, std::vector& user_solution); +template void uncrush_dual_solution(const user_problem_t& user_problem, + const lp_problem_t& problem, + const std::vector& y, + const std::vector& z, + std::vector& user_y, + std::vector& user_z); + template void uncrush_solution(const presolve_info_t& presolve_info, const std::vector& crushed_x, const std::vector& crushed_z, diff --git a/cpp/src/dual_simplex/presolve.hpp b/cpp/src/dual_simplex/presolve.hpp index e6eb542ac..947c637cb 100644 --- a/cpp/src/dual_simplex/presolve.hpp +++ b/cpp/src/dual_simplex/presolve.hpp @@ -115,6 +115,14 @@ void uncrush_primal_solution(const user_problem_t& user_problem, const std::vector& solution, std::vector& user_solution); +template +void uncrush_dual_solution(const user_problem_t& user_problem, + const lp_problem_t& problem, + const std::vector& y, + const std::vector& z, + std::vector& user_y, + std::vector& user_z); + template void uncrush_solution(const presolve_info_t& presolve_info, const std::vector& crushed_x, diff --git a/cpp/src/dual_simplex/scaling.cpp b/cpp/src/dual_simplex/scaling.cpp index 3a19ec613..50f322ee5 100644 --- a/cpp/src/dual_simplex/scaling.cpp +++ b/cpp/src/dual_simplex/scaling.cpp @@ -88,7 +88,7 @@ void unscale_solution(const std::vector& column_scaling, unscaled_z.resize(n); for (i_t j = 0; j < n; ++j) { unscaled_x[j] = scaled_x[j] / column_scaling[j]; - unscaled_z[j] = scaled_z[j] / column_scaling[j]; + unscaled_z[j] = scaled_z[j] * column_scaling[j]; } } diff --git a/cpp/src/dual_simplex/solve.cpp b/cpp/src/dual_simplex/solve.cpp index 29aecb21f..e1497fae6 100644 --- a/cpp/src/dual_simplex/solve.cpp +++ b/cpp/src/dual_simplex/solve.cpp @@ -252,8 +252,8 @@ lp_status_t solve_linear_program(const user_problem_t& user_problem, lp_status_t status = solve_linear_program_advanced( original_lp, start_time, settings, lp_solution, vstatus, edge_norms); uncrush_primal_solution(user_problem, original_lp, lp_solution.x, solution.x); - uncrush_primal_solution(user_problem, original_lp, lp_solution.z, solution.z); - solution.y = lp_solution.y; + uncrush_dual_solution( + user_problem, original_lp, lp_solution.y, lp_solution.z, solution.y, solution.z); solution.objective = lp_solution.objective; solution.user_objective = lp_solution.user_objective; solution.iterations = lp_solution.iterations; diff --git a/cpp/tests/dual_simplex/unit_tests/solve.cpp b/cpp/tests/dual_simplex/unit_tests/solve.cpp index 3fbfc3766..0743c08d7 100644 --- a/cpp/tests/dual_simplex/unit_tests/solve.cpp +++ b/cpp/tests/dual_simplex/unit_tests/solve.cpp @@ -262,4 +262,74 @@ TEST(dual_simplex, empty_columns) EXPECT_NEAR(solution.x[8], 0, 1e-6); } +TEST(dual_simplex, dual_variable_greater_than) +{ + // minimize 3*x0 + 2 * x1 + // subject to x0 + x1 >= 1 + // x0 + 2x1 >= 3 + // x0, x1 >= 0 + + cuopt::linear_programming::dual_simplex::user_problem_t user_problem; + constexpr int m = 2; + constexpr int n = 2; + constexpr int nz = 4; + + user_problem.num_rows = m; + user_problem.num_cols = n; + user_problem.objective.resize(n); + user_problem.objective[0] = 3.0; + user_problem.objective[1] = 2.0; + user_problem.A.m = m; + user_problem.A.n = n; + user_problem.A.nz_max = nz; + user_problem.A.reallocate(nz); + user_problem.A.col_start.resize(n + 1); + user_problem.A.col_start[0] = 0; // x0 start + user_problem.A.col_start[1] = 2; + user_problem.A.col_start[2] = 4; + + int nnz = 0; + user_problem.A.i[nnz] = 0; + user_problem.A.x[nnz++] = 1.0; + user_problem.A.i[nnz] = 1; + user_problem.A.x[nnz++] = 1.0; + user_problem.A.i[nnz] = 0; + user_problem.A.x[nnz++] = 1.0; + user_problem.A.i[nnz] = 1; + user_problem.A.x[nnz++] = 2.0; + user_problem.A.print_matrix(); + EXPECT_EQ(nnz, nz); + + user_problem.rhs.resize(m); + user_problem.rhs[0] = 1.0; + user_problem.rhs[1] = 3.0; + + user_problem.row_sense.resize(m); + user_problem.row_sense[0] = 'G'; + user_problem.row_sense[1] = 'G'; + + user_problem.lower.resize(n); + user_problem.lower[0] = 0.0; + user_problem.lower[1] = 0.0; + + user_problem.upper.resize(n); + user_problem.upper[0] = dual_simplex::inf; + user_problem.upper[1] = dual_simplex::inf; + + user_problem.num_range_rows = 0; + user_problem.problem_name = "dual_variable_greater_than"; + + dual_simplex::simplex_solver_settings_t settings; + dual_simplex::lp_solution_t solution(user_problem.num_rows, user_problem.num_cols); + EXPECT_EQ((dual_simplex::solve_linear_program(user_problem, settings, solution)), + dual_simplex::lp_status_t::OPTIMAL); + EXPECT_NEAR(solution.objective, 3.0, 1e-6); + EXPECT_NEAR(solution.x[0], 0.0, 1e-6); + EXPECT_NEAR(solution.x[1], 1.5, 1e-6); + EXPECT_NEAR(solution.y[0], 0.0, 1e-6); + EXPECT_NEAR(solution.y[1], 1.0, 1e-6); + EXPECT_NEAR(solution.z[0], 2.0, 1e-6); + EXPECT_NEAR(solution.z[1], 0.0, 1e-6); +} + } // namespace cuopt::linear_programming::dual_simplex::test