diff --git a/cpp/cuopt_cli.cpp b/cpp/cuopt_cli.cpp index 4f6422dfc..a9aeaef77 100644 --- a/cpp/cuopt_cli.cpp +++ b/cpp/cuopt_cli.cpp @@ -104,7 +104,7 @@ int run_single_file(const std::string& file_path, cuopt::mps_parser::mps_data_model_t mps_data_model; bool parsing_failed = false; { - CUOPT_LOG_INFO("Running file %s", base_filename.c_str()); + CUOPT_LOG_INFO("Reading file %s", base_filename.c_str()); try { mps_data_model = cuopt::mps_parser::parse_mps(file_path, input_mps_strict); } catch (const std::logic_error& e) { diff --git a/cpp/include/cuopt/linear_programming/optimization_problem.hpp b/cpp/include/cuopt/linear_programming/optimization_problem.hpp index 6d9ba5cda..6f8f95973 100644 --- a/cpp/include/cuopt/linear_programming/optimization_problem.hpp +++ b/cpp/include/cuopt/linear_programming/optimization_problem.hpp @@ -311,9 +311,13 @@ class optimization_problem_t { */ void write_to_mps(const std::string& mps_file_path); + /* Print scaling information */ + void print_scaling_information() const; + i_t get_n_variables() const; i_t get_n_constraints() const; i_t get_nnz() const; + i_t get_n_integers() const; raft::handle_t const* get_handle_ptr() const noexcept; const rmm::device_uvector& get_constraint_matrix_values() const; rmm::device_uvector& get_constraint_matrix_values(); diff --git a/cpp/src/linear_programming/optimization_problem.cu b/cpp/src/linear_programming/optimization_problem.cu index 5847c503b..bcd2f9c4a 100644 --- a/cpp/src/linear_programming/optimization_problem.cu +++ b/cpp/src/linear_programming/optimization_problem.cu @@ -16,6 +16,7 @@ */ #include +#include #include #include @@ -273,6 +274,20 @@ i_t optimization_problem_t::get_nnz() const return A_.size(); } +template +i_t optimization_problem_t::get_n_integers() const +{ + i_t n_integers = 0; + if (get_n_variables() != 0) { + auto enum_variable_types = cuopt::host_copy(get_variable_types()); + + for (size_t i = 0; i < enum_variable_types.size(); ++i) { + if (enum_variable_types[i] == var_t::INTEGER) { n_integers++; } + } + } + return n_integers; +} + template raft::handle_t const* optimization_problem_t::get_handle_ptr() const noexcept { @@ -583,6 +598,84 @@ void optimization_problem_t::write_to_mps(const std::string& mps_file_ cuopt::mps_parser::write_mps(data_model_view, mps_file_path); } +template +void optimization_problem_t::print_scaling_information() const +{ + std::vector constraint_matrix_values = cuopt::host_copy(get_constraint_matrix_values()); + std::vector constraint_rhs = cuopt::host_copy(get_constraint_bounds()); + std::vector objective_coefficients = cuopt::host_copy(get_objective_coefficients()); + std::vector variable_lower_bounds = cuopt::host_copy(get_variable_lower_bounds()); + std::vector variable_upper_bounds = cuopt::host_copy(get_variable_upper_bounds()); + std::vector constraint_lower_bounds = cuopt::host_copy(get_constraint_lower_bounds()); + std::vector constraint_upper_bounds = cuopt::host_copy(get_constraint_upper_bounds()); + + auto findMaxAbs = [](const std::vector& vec) -> f_t { + if (vec.empty()) { return 0.0; } + const f_t inf = std::numeric_limits::infinity(); + + const size_t sz = vec.size(); + f_t max_abs_val = 0.0; + for (size_t i = 0; i < sz; ++i) { + const f_t val = std::abs(vec[i]); + if (val < inf) { max_abs_val = std::max(max_abs_val, val); } + } + return max_abs_val; + }; + + auto findMinAbs = [](const std::vector& vec) -> f_t { + if (vec.empty()) { return 0.0; } + const size_t sz = vec.size(); + const f_t inf = std::numeric_limits::infinity(); + f_t min_abs_val = inf; + for (size_t i = 0; i < sz; ++i) { + const f_t val = std::abs(vec[i]); + if (val > 0.0) { min_abs_val = std::min(min_abs_val, val); } + } + return min_abs_val < inf ? min_abs_val : 0.0; + }; + + f_t A_max = findMaxAbs(constraint_matrix_values); + f_t A_min = findMinAbs(constraint_matrix_values); + f_t b_max = findMaxAbs(constraint_rhs); + f_t b_min = findMinAbs(constraint_rhs); + f_t c_max = findMaxAbs(objective_coefficients); + f_t c_min = findMinAbs(objective_coefficients); + f_t x_lower_max = findMaxAbs(variable_lower_bounds); + f_t x_lower_min = findMinAbs(variable_lower_bounds); + f_t x_upper_max = findMaxAbs(variable_upper_bounds); + f_t x_upper_min = findMinAbs(variable_upper_bounds); + f_t cstr_lower_max = findMaxAbs(constraint_lower_bounds); + f_t cstr_lower_min = findMinAbs(constraint_lower_bounds); + f_t cstr_upper_max = findMaxAbs(constraint_upper_bounds); + f_t cstr_upper_min = findMinAbs(constraint_upper_bounds); + + f_t rhs_max = std::max(b_max, std::max(cstr_lower_max, cstr_upper_max)); + f_t rhs_min = std::min(b_min, std::min(cstr_lower_min, cstr_upper_min)); + + f_t bound_max = std::max(x_upper_max, x_lower_max); + f_t bound_min = std::min(x_upper_min, x_lower_min); + + CUOPT_LOG_INFO("Problem scaling:"); + CUOPT_LOG_INFO("Objective coefficents range: [%.0e, %.0e]", c_min, c_max); + CUOPT_LOG_INFO("Constraint matrix coefficients range: [%.0e, %.0e]", A_min, A_max); + CUOPT_LOG_INFO("Constraint rhs / bounds range: [%.0e, %.0e]", rhs_min, rhs_max); + CUOPT_LOG_INFO("Variable bounds range: [%.0e, %.0e]", bound_min, bound_max); + + auto safelog10 = [](f_t x) { return x > 0 ? std::log10(x) : 0.0; }; + + f_t obj_range = safelog10(c_max) - safelog10(c_min); + f_t A_range = safelog10(A_max) - safelog10(A_min); + f_t rhs_range = safelog10(rhs_max) - safelog10(rhs_min); + f_t bound_range = safelog10(bound_max) - safelog10(bound_min); + + if (obj_range >= 6.0 || A_range >= 6.0 || rhs_range >= 6.0 || bound_range >= 6.0) { + CUOPT_LOG_INFO( + "Warning: input problem contains a large range of coefficients: consider reformulating to " + "avoid numerical difficulties."); + } + CUOPT_LOG_INFO(""); +} + // NOTE: Explicitly instantiate all types here in order to avoid linker error #if MIP_INSTANTIATE_FLOAT template class optimization_problem_t; diff --git a/cpp/src/linear_programming/solve.cu b/cpp/src/linear_programming/solve.cu index 95c7f3333..f8f88f8f8 100644 --- a/cpp/src/linear_programming/solve.cu +++ b/cpp/src/linear_programming/solve.cu @@ -813,6 +813,14 @@ optimization_problem_solution_t solve_lp(optimization_problem_t::check_initial_solution_representation(op_problem, settings); } + CUOPT_LOG_INFO( + "Solving a problem with %d constraints, %d variables (%d integers), and %d nonzeros", + op_problem.get_n_constraints(), + op_problem.get_n_variables(), + 0, + op_problem.get_nnz()); + op_problem.print_scaling_information(); + // Check for crossing bounds. Return infeasible if there are any if (problem_checking_t::has_crossing_bounds(op_problem)) { return optimization_problem_solution_t(pdlp_termination_status_t::PrimalInfeasible, @@ -851,12 +859,6 @@ optimization_problem_solution_t solve_lp(optimization_problem_t, bool> third_party_presolve_t papilo_problem = build_papilo_problem(op_problem, category); - CUOPT_LOG_INFO("Unpresolved problem: %d constraints, %d variables, %d nonzeros", + CUOPT_LOG_INFO("Original problem: %d constraints, %d variables, %d nonzeros", papilo_problem.getNRows(), papilo_problem.getNCols(), papilo_problem.getConstraintMatrix().getNnz()); diff --git a/cpp/src/mip/solve.cu b/cpp/src/mip/solve.cu index 547abee49..974cde0bc 100644 --- a/cpp/src/mip/solve.cu +++ b/cpp/src/mip/solve.cu @@ -96,11 +96,7 @@ mip_solution_t run_mip(detail::problem_t& problem, } // problem contains unpreprocessed data detail::problem_t scaled_problem(problem); - CUOPT_LOG_INFO("Solving a problem with %d constraints %d variables (%d integers) and %d nonzeros", - problem.n_constraints, - problem.n_variables, - problem.n_integer_vars, - problem.nnz); + CUOPT_LOG_INFO("Objective offset %f scaling_factor %f", problem.presolve_data.objective_offset, problem.presolve_data.objective_scaling_factor); @@ -181,6 +177,14 @@ mip_solution_t solve_mip(optimization_problem_t& op_problem, problem_checking_t::check_problem_representation(op_problem); problem_checking_t::check_initial_solution_representation(op_problem, settings); + CUOPT_LOG_INFO( + "Solving a problem with %d constraints, %d variables (%d integers), and %d nonzeros", + op_problem.get_n_constraints(), + op_problem.get_n_variables(), + op_problem.get_n_integers(), + op_problem.get_nnz()); + op_problem.print_scaling_information(); + // Check for crossing bounds. Return infeasible if there are any if (problem_checking_t::has_crossing_bounds(op_problem)) { return mip_solution_t(mip_termination_status_t::Infeasible,