diff --git a/cpp/src/barrier/translate_soc.hpp b/cpp/src/barrier/translate_soc.hpp index 3209dee3c6..1da776188f 100644 --- a/cpp/src/barrier/translate_soc.hpp +++ b/cpp/src/barrier/translate_soc.hpp @@ -659,9 +659,15 @@ void convert_quadratic_constraints_to_second_order_cones( "Quadratic constraint '%s' is non-convex (Q matrix is indefinite)", qc.constraint_row_name.c_str()); - // Since q_nnz >= 1 is enforced above, Q is nonzero and rank must be >= 1. - // (A nonzero Q entry produces a nonzero H diagonal or off-diagonal, guaranteeing rank > 0.) - assert(rank >= 1); + // q_nnz >= 1 implies H is nonzero, but diagonal LDLT may still return rank 0 (e.g. cross-only + // indefinite H with zero diagonals). Reject before building a degenerate r=0 SOC lift. + cuopt_expects(rank >= 1, + error_type_t::ValidationError, + "Quadratic constraint '%s' is non-convex or could not be converted to a " + "second-order cone (LDLT rank %d; Q may be indefinite or have zero diagonal " + "with cross terms)", + qc.constraint_row_name.c_str(), + static_cast(rank)); // Step 4: Build standard SOC of dimension rank + 2. // New variables: y_0,...,y_{r-1}, s_0 (head), s_{r+1} (tail) diff --git a/cpp/src/dual_simplex/right_looking_lu.cpp b/cpp/src/dual_simplex/right_looking_lu.cpp index 00f2fefb62..f387cadc53 100644 --- a/cpp/src/dual_simplex/right_looking_lu.cpp +++ b/cpp/src/dual_simplex/right_looking_lu.cpp @@ -13,7 +13,6 @@ #include #include -#include using cuopt::ins_vector; @@ -1684,7 +1683,8 @@ i_t right_looking_ldlt(const csc_matrix_t& A, f_t& work_estimate) { raft::common::nvtx::range scope("LU::right_looking_ldlt"); - const i_t n = A.n; + const i_t n = A.n; + const i_t input_nnz = A.nnz(); assert(A.m == n); symmetric_trailing_matrix_t trailing_matrix(A); @@ -1720,7 +1720,12 @@ i_t right_looking_ldlt(const csc_matrix_t& A, f_t pivot_val = 0; trailing_matrix.symmetric_markowitz_search(pivot_tol, pivot_p, pivot_val); - if (pivot_p == -1) { break; } // No acceptable pivot found (remaining diagonals are zero/tiny) + if (pivot_p == -1) { + // No acceptable diagonal pivot. Zero matrix is rank 0; nonzero cross-only indefinite + // matrices (e.g. [[0, a]; [a, 0]]) also stall here and must not be treated as PSD. + if (pivots == 0 && input_nnz > 0) { return INDEFINITE_MATRIX_RETURN; } + break; + } // Check for indefiniteness: a negative pivot means the matrix is not PSD if (pivot_val < 0) { return INDEFINITE_MATRIX_RETURN; } diff --git a/cpp/src/dual_simplex/right_looking_lu.hpp b/cpp/src/dual_simplex/right_looking_lu.hpp index 9b2814dfde..84b2604b53 100644 --- a/cpp/src/dual_simplex/right_looking_lu.hpp +++ b/cpp/src/dual_simplex/right_looking_lu.hpp @@ -43,7 +43,8 @@ i_t right_looking_lu_row_permutation_only(const csc_matrix_t& A, // D = diagonal factor (length = rank) // Returns: // rank >= 0: number of successful pivots with D(k,k) >= pivot_tol (PSD case). -// INDEFINITE_MATRIX_RETURN (-4): a negative pivot was encountered (matrix is not PSD). +// INDEFINITE_MATRIX_RETURN (-4): a negative pivot was encountered, or the matrix is nonzero +// but no acceptable diagonal pivot was found (matrix is not PSD). // CONCURRENT_HALT_RETURN (-2): concurrent halt requested. // TIME_LIMIT_RETURN (-3): time limit exceeded. template diff --git a/cpp/src/pdlp/solve.cu b/cpp/src/pdlp/solve.cu index 5d29e97fa0..3ea3d42dc2 100644 --- a/cpp/src/pdlp/solve.cu +++ b/cpp/src/pdlp/solve.cu @@ -2124,84 +2124,7 @@ cuopt::linear_programming::optimization_problem_t mps_data_model_to_op error_type_t::ValidationError, "handle_ptr must not be null for GPU-backed problem construction"); cuopt::linear_programming::optimization_problem_t op_problem(handle_ptr); - op_problem.set_maximize(data_model.get_sense()); - - if (data_model.get_constraint_matrix_values().size() != 0) { - op_problem.set_csr_constraint_matrix(data_model.get_constraint_matrix_values().data(), - data_model.get_constraint_matrix_values().size(), - data_model.get_constraint_matrix_indices().data(), - data_model.get_constraint_matrix_indices().size(), - data_model.get_constraint_matrix_offsets().data(), - data_model.get_constraint_matrix_offsets().size()); - } else { - // Set empty constraint matrix - std::vector offsets(1, 0); - op_problem.set_csr_constraint_matrix(nullptr, 0, nullptr, 0, offsets.data(), 1); - } - - if (data_model.get_constraint_bounds().size() != 0) { - op_problem.set_constraint_bounds(data_model.get_constraint_bounds().data(), - data_model.get_constraint_bounds().size()); - } - if (data_model.get_objective_coefficients().size() != 0) { - op_problem.set_objective_coefficients(data_model.get_objective_coefficients().data(), - data_model.get_objective_coefficients().size()); - } - op_problem.set_objective_scaling_factor(data_model.get_objective_scaling_factor()); - op_problem.set_objective_offset(data_model.get_objective_offset()); - if (data_model.get_variable_lower_bounds().size() != 0) { - op_problem.set_variable_lower_bounds(data_model.get_variable_lower_bounds().data(), - data_model.get_variable_lower_bounds().size()); - } - if (data_model.get_variable_upper_bounds().size() != 0) { - op_problem.set_variable_upper_bounds(data_model.get_variable_upper_bounds().data(), - data_model.get_variable_upper_bounds().size()); - } - if (data_model.get_variable_types().size() != 0) { - std::vector enum_variable_types(data_model.get_variable_types().size()); - std::transform(data_model.get_variable_types().cbegin(), - data_model.get_variable_types().cend(), - enum_variable_types.begin(), - detail::char_to_var_type); - op_problem.set_variable_types(enum_variable_types.data(), enum_variable_types.size()); - } - - if (data_model.get_row_types().size() != 0) { - op_problem.set_row_types(data_model.get_row_types().data(), data_model.get_row_types().size()); - } - if (data_model.get_constraint_lower_bounds().size() != 0) { - op_problem.set_constraint_lower_bounds(data_model.get_constraint_lower_bounds().data(), - data_model.get_constraint_lower_bounds().size()); - } - if (data_model.get_constraint_upper_bounds().size() != 0) { - op_problem.set_constraint_upper_bounds(data_model.get_constraint_upper_bounds().data(), - data_model.get_constraint_upper_bounds().size()); - } - - if (data_model.get_objective_name().size() != 0) { - op_problem.set_objective_name(data_model.get_objective_name()); - } - auto problem_name = data_model.get_problem_name(); - op_problem.set_problem_name(problem_name); - if (data_model.get_variable_names().size() != 0) { - op_problem.set_variable_names(data_model.get_variable_names()); - } - if (data_model.get_row_names().size() != 0) { - op_problem.set_row_names(data_model.get_row_names()); - } - - if (data_model.get_quadratic_objective_values().size() != 0) { - const std::vector Q_values = data_model.get_quadratic_objective_values(); - const std::vector Q_indices = data_model.get_quadratic_objective_indices(); - const std::vector Q_offsets = data_model.get_quadratic_objective_offsets(); - op_problem.set_quadratic_objective_matrix(Q_values.data(), - Q_values.size(), - Q_indices.data(), - Q_indices.size(), - Q_offsets.data(), - Q_offsets.size()); - } - + populate_from_mps_data_model(&op_problem, data_model); return op_problem; } diff --git a/cpp/src/pdlp/translate.hpp b/cpp/src/pdlp/translate.hpp index be832e082d..4300e9c13e 100644 --- a/cpp/src/pdlp/translate.hpp +++ b/cpp/src/pdlp/translate.hpp @@ -44,8 +44,10 @@ static dual_simplex::user_problem_t cuopt_problem_to_user_problem( csr_A.x = std::move(A_values); csr_A.j = std::move(A_indices); csr_A.row_start = std::move(A_offsets); - - csr_A.to_compressed_col(user_problem.A); + if (m == 0) { + csr_A.row_start.resize(1); + csr_A.row_start[0] = 0; + } user_problem.rhs.resize(m); user_problem.row_sense.resize(m); @@ -96,6 +98,13 @@ static dual_simplex::user_problem_t cuopt_problem_to_user_problem( user_problem.Q_indices = problem.get_quadratic_objective_indices(); user_problem.Q_values = problem.get_quadratic_objective_values(); + if (problem.has_quadratic_constraints()) { + detail::convert_quadratic_constraints_to_second_order_cones( + n, problem.get_quadratic_constraints(), csr_A, user_problem); + } + + csr_A.to_compressed_col(user_problem.A); + return user_problem; } diff --git a/cpp/tests/dual_simplex/unit_tests/right_looking_ldlt.cpp b/cpp/tests/dual_simplex/unit_tests/right_looking_ldlt.cpp index 7c91a1f72f..34f2c706d5 100644 --- a/cpp/tests/dual_simplex/unit_tests/right_looking_ldlt.cpp +++ b/cpp/tests/dual_simplex/unit_tests/right_looking_ldlt.cpp @@ -512,6 +512,26 @@ TEST(right_looking_ldlt, indefinite_2x2) EXPECT_EQ(result, INDEFINITE_MATRIX_RETURN); } +// Test 14b: Cross-only indefinite matrix (zero diagonals, no 1x1 pivot). +// A = [0, 2; 2, 0] has eigenvalues +2 and -2 (indefinite). +TEST(right_looking_ldlt, indefinite_cross_only_2x2) +{ + const int n = 2; + std::vector dense = {0.0, 2.0, 2.0, 0.0}; + auto A = dense_to_lower_csc(n, dense); + + simplex_solver_settings_t settings; + std::vector perm; + csc_matrix_t L(n, n, 1); + std::vector D; + double work_estimate = 0; + double start_time = tic(); + + int result = right_looking_ldlt(A, settings, 1e-12, start_time, perm, L, D, work_estimate); + + EXPECT_EQ(result, INDEFINITE_MATRIX_RETURN); +} + // Test 15: Larger indefinite matrix (4x4 with mixed eigenvalues). // A = [2, 1, 0, 1; 1, 0, 1, 0; 0, 1, 2, -1; 1, 0, -1, -1] TEST(right_looking_ldlt, indefinite_4x4) diff --git a/cpp/tests/linear_programming/unit_tests/solution_interface_test.cu b/cpp/tests/linear_programming/unit_tests/solution_interface_test.cu index f7f164900a..7ce05e60d6 100644 --- a/cpp/tests/linear_programming/unit_tests/solution_interface_test.cu +++ b/cpp/tests/linear_programming/unit_tests/solution_interface_test.cu @@ -400,6 +400,30 @@ TEST_F(SolutionInterfaceTest, mps_data_model_to_optimization_problem) } } +TEST_F(SolutionInterfaceTest, mps_data_model_to_optimization_problem_quadratic_constraints) +{ + auto mps_data = cuopt::linear_programming::io::read_lp_from_string(R"LP( +Minimize + obj: x + y +Subject To + q0: [ 4 x * y ] <= 0.5 +Bounds + -1 <= x <= 1 + -1 <= y <= 1 +End +)LP"); + ASSERT_TRUE(mps_data.has_quadratic_constraints()); + ASSERT_EQ(mps_data.get_quadratic_constraints().size(), 1u); + + raft::handle_t handle; + auto problem = mps_data_model_to_optimization_problem(&handle, mps_data); + + ASSERT_TRUE(problem.has_quadratic_constraints()); + ASSERT_EQ(problem.get_quadratic_constraints().size(), 1u); + EXPECT_EQ(problem.get_quadratic_constraints()[0].constraint_row_name, "q0"); + EXPECT_NEAR(problem.get_quadratic_constraints()[0].rhs_value, 0.5, 1e-9); +} + // ============================================================================= // Solution conversion tests (hand-constructed, known values) // ============================================================================= diff --git a/cpp/tests/socp/general_quadratic_test.cu b/cpp/tests/socp/general_quadratic_test.cu index 2918d80a7d..935f25359f 100644 --- a/cpp/tests/socp/general_quadratic_test.cu +++ b/cpp/tests/socp/general_quadratic_test.cu @@ -5,13 +5,20 @@ */ /* clang-format on */ +#include #include #include +#include +#include #include +#include +#include +#include #include #include #include +#include #include #include @@ -207,6 +214,44 @@ TEST(general_quadratic, rejects_non_convex) cuopt::logic_error); } +// End-to-end: cross-only indefinite Q (issue #1434). H = [[0, 2]; [2, 0]] has zero diagonals so +// LDLT returns rank 0 without a negative pivot; must still be rejected as non-convex via +// solve_qcqp. +TEST(general_quadratic, rejects_cross_only_indefinite) +{ + raft::handle_t handle{}; + init_handler(&handle); + + // H = [[0, 2]; [2, 0]] from [ 4 x * y ] in LP bracket notation. + auto problem = cuopt::linear_programming::io::read_lp_from_string(R"LP( +Minimize + obj: x + y +Subject To + q0: [ 4 x * y ] <= 0.5 +Bounds + -1 <= x <= 1 + -1 <= y <= 1 +End +)LP"); + + ASSERT_TRUE(problem.has_quadratic_constraints()); + ASSERT_EQ(problem.get_quadratic_constraints().size(), 1u); + + cuopt::linear_programming::optimization_problem_t op_problem(&handle); + cuopt::linear_programming::populate_from_mps_data_model(&op_problem, problem); + ASSERT_TRUE(op_problem.has_quadratic_constraints()) + << "populate_from_mps_data_model dropped quadratic constraints"; + ASSERT_EQ(op_problem.get_quadratic_constraints().size(), 1u); + + cuopt::linear_programming::pdlp_solver_settings_t settings; + auto solution = cuopt::linear_programming::solve_lp(op_problem, settings); + + const auto error_status = solution.get_error_status(); + EXPECT_EQ(error_status.get_error_type(), cuopt::error_type_t::ValidationError); + EXPECT_THAT(error_status.what(), testing::HasSubstr("non-convex")); + EXPECT_THAT(error_status.what(), testing::HasSubstr("q0")); +} + // Test: rank-deficient PSD Q (e.g., Q = v*v^T with v = [1, 1]) // minimize x0 + x1 // subject to (x0 + x1)^2 <= 4 (i.e., |x0 + x1| <= 2)