Add support for quadratic programming to the barrier method (#500)

This PR adds support for quadratic programs (QPs) using the barrier method. Only quadratic terms in the objective (pure QPs) are supported.  Quadratic constraints (or QCQPs) are not supported at this time. 

The following changes are included in this PR: 
- The barrier method is modified to support quadratic objectives. LPs remain unchanged.
- When Q matrix is diagonal,  normal equations are formed, otherwise augmented system is formed
- C and C++ APIs for setting up quadratic objective 
- Python API for setting Q matrix in CSR format
- MPS and QPS/SIF formats are supported and can be run directly with cuopt_cli

Authors:
  - Chris Maes (https://github.com/chris-maes)
  - Ramakrishnap (https://github.com/rgsl888prabhu)
  - Rajesh Gandham (https://github.com/rg20)

Approvers:
  - Rajesh Gandham (https://github.com/rg20)
  - Trevor McKay (https://github.com/tmckayus)

URL: #500

Author: chris-maes

benchmarks/linear_programming/run_mps_files.sh

@@ -127,6 +127,11 @@ while [[ $# -gt 0 ]]; do
             TIME_LIMIT="$2"
             shift 2
             ;;
+        --iteration-limit)
+            echo "ITERATION_LIMIT: $2"
+            ITERATION_LIMIT="$2"
+            shift 2
+            ;;
         --output-dir)
             echo "OUTPUT_DIR: $2"
             OUTPUT_DIR="$2"
@@ -152,6 +157,11 @@ while [[ $# -gt 0 ]]; do
             NUM_CPU_THREADS="$2"
             shift 2
             ;;
+        --method)
+            echo "METHOD: $2"
+            METHOD="$2"
+            shift 2
+            ;;
         --presolve)
             echo "PRESOLVE: $2"
             PRESOLVE="$2"
@@ -200,7 +210,6 @@ RELAXATION=${RELAXATION:-false}
 MIP_HEURISTICS_ONLY=${MIP_HEURISTICS_ONLY:-false}
 WRITE_LOG_FILE=${WRITE_LOG_FILE:-false}
 NUM_CPU_THREADS=${NUM_CPU_THREADS:-1}
-PRESOLVE=${PRESOLVE:-true}
 BATCH_NUM=${BATCH_NUM:-0}
 N_BATCHES=${N_BATCHES:-1}
 LOG_TO_CONSOLE=${LOG_TO_CONSOLE:-true}
@@ -281,7 +290,8 @@ if [[ -n "$MODEL_LIST" && -f "$MODEL_LIST" ]]; then
         # Skip empty lines
         [[ -z "$line" ]] && continue
         # If the line does not end with .mps, append it
-        if [[ "$line" != *.mps ]]; then
+        # ignore if it ends with .SIF
+        if [[ "$line" != *.mps && "$line" != *.SIF ]]; then
             echo "${line}.mps" >> "$TMP_MODEL_LIST"
         else
             echo "$line" >> "$TMP_MODEL_LIST"
@@ -314,7 +324,10 @@ if [[ -n "$MODEL_LIST" ]]; then
         exit 1
     fi
 else
-    mapfile -t mps_files < <(ls "$MPS_DIR"/*.mps)
+    # Gather both .mps and .SIF files in the directory
+    mapfile -t mps_files < <(ls "$MPS_DIR"/*.mps "$MPS_DIR"/*.SIF 2>/dev/null)
+
+    echo "Found ${#mps_files[@]} .mps and .SIF files in $MPS_DIR"
 fi
 
 # Calculate batch size and start/end indices
@@ -375,6 +388,9 @@ worker() {
 
         # Build arguments string
         args=""
+        if [ -n "$ITERATION_LIMIT" ]; then
+            args="$args --iteration-limit $ITERATION_LIMIT"
+        fi
         if [ -n "$NUM_CPU_THREADS" ]; then
             args="$args --num-cpu-threads $NUM_CPU_THREADS"
         fi
@@ -391,7 +407,12 @@ worker() {
             args="$args --num-gpus $GPUS_PER_INSTANCE"
         fi
         args="$args --log-to-console $LOG_TO_CONSOLE"
-        args="$args --presolve $PRESOLVE"
+        if [ -n "$PRESOLVE" ]; then
+            args="$args --presolve $PRESOLVE"
+        fi
+        if [ -n "$METHOD" ]; then
+            args="$args --method $METHOD"
+        fi
 
         CUDA_VISIBLE_DEVICES=$gpu_devices cuopt_cli "$mps_file" --time-limit $TIME_LIMIT $args
     done

cpp/include/cuopt/linear_programming/cuopt_c.h

@@ -114,9 +114,9 @@ cuopt_int_t cuOptReadProblem(const char* filename, cuOptOptimizationProblem* pro
 /** @brief Create an optimization problem of the form
  *
  * @verbatim
- *                minimize/maximize  cᵀx + offset
- *                  subject to       A x {=, ≤, ≥} b
- *                                   l ≤ x ≤ u
+ *                minimize/maximize  c^T x + offset
+ *                  subject to       A x {=, <=, >=} b
+ *                                   l <= x <= u
  *                                   x_i integer for some i
  * @endverbatim
  *
@@ -238,6 +238,160 @@ cuopt_int_t cuOptCreateRangedProblem(cuopt_int_t num_constraints,
                                      const char* variable_types,
                                      cuOptOptimizationProblem* problem_ptr);
 
+/** @brief Create an optimization problem of the form
+ *
+ * @verbatim
+ *                minimize/maximize  c^T x + x^T Q x + offset
+ *                  subject to       A x {=, <=, >=} b
+ *                                   l ≤ x ≤ u
+ * @endverbatim
+ *
+ * @param[in] num_constraints The number of constraints
+ * @param[in] num_variables The number of variables
+ * @param[in] objective_sense The objective sense (CUOPT_MINIMIZE for
+ *            minimization or CUOPT_MAXIMIZE for maximization)
+ * @param[in] objective_offset An offset to add to the linear objective
+ * @param[in] objective_coefficients A pointer to an array of type cuopt_float_t
+ *            of size num_variables containing the coefficients of the linear objective
+ * @param[in] quadratic_objective_matrix_row_offsets A pointer to an array of type
+ *            cuopt_int_t of size num_variables + 1. quadratic_objective_matrix_row_offsets[i] is
+ * the index of the first non-zero element of the i-th row of the quadratic objective matrix in
+ *            quadratic_objective_matrix_column_indices and
+ * quadratic_objective_matrix_coefficent_values. This is part of the compressed sparse row
+ * representation of the quadratic objective matrix.
+ * @param[in] quadratic_objective_matrix_column_indices A pointer to an array of type
+ *            cuopt_int_t of size quadratic_objective_matrix_row_offsets[num_variables] containing
+ *            the column indices of the non-zero elements of the quadratic objective matrix.
+ *            This is part of the compressed sparse row representation of the quadratic objective
+ * matrix.
+ * @param[in] quadratic_objective_matrix_coefficent_values A pointer to an array of type
+ *            cuopt_float_t of size quadratic_objective_matrix_row_offsets[num_variables] containing
+ *            the values of the non-zero elements of the quadratic objective matrix.
+ * @param[in] constraint_matrix_row_offsets A pointer to an array of type
+ *            cuopt_int_t of size num_constraints + 1. constraint_matrix_row_offsets[i] is the
+ *            index of the first non-zero element of the i-th constraint in
+ *            constraint_matrix_column_indices and constraint_matrix_coefficent_values. This is
+ *            part of the compressed sparse row representation of the constraint matrix
+ * @param[in] constraint_matrix_column_indices A pointer to an array of type
+ *            cuopt_int_t of size constraint_matrix_row_offsets[num_constraints] containing
+ *            the column indices of the non-zero elements of the constraint matrix. This is
+ *            part of the compressed sparse row representation of the constraint matrix
+ * @param[in] constraint_matrix_coefficent_values A pointer to an array of type
+ *            cuopt_float_t of size constraint_matrix_row_offsets[num_constraints] containing
+ *            the values of the non-zero elements of the constraint matrix. This is
+ *            part of the compressed sparse row representation of the constraint matrix
+ * @param[in] constraint_sense A pointer to an array of type char of size
+ *            num_constraints containing the sense of the constraints (CUOPT_LESS_THAN,
+ *            CUOPT_GREATER_THAN, or CUOPT_EQUAL)
+ * @param[in] rhs A pointer to an array of type cuopt_float_t of size num_constraints
+ *            containing the right-hand side of the constraints
+ * @param[in] lower_bounds A pointer to an array of type cuopt_float_t of size num_variables
+ *            containing the lower bounds of the variables
+ * @param[in] upper_bounds A pointer to an array of type cuopt_float_t of size num_variables
+ *            containing the upper bounds of the variables
+ * @param[out] problem_ptr Pointer to store the created optimization problem
+ * @return CUOPT_SUCCESS if successful, CUOPT_ERROR otherwise
+ */
+cuopt_int_t cuOptCreateQuadraticProblem(
+  cuopt_int_t num_constraints,
+  cuopt_int_t num_variables,
+  cuopt_int_t objective_sense,
+  cuopt_float_t objective_offset,
+  const cuopt_float_t* objective_coefficients,
+  const cuopt_int_t* quadratic_objective_matrix_row_offsets,
+  const cuopt_int_t* quadratic_objective_matrix_column_indices,
+  const cuopt_float_t* quadratic_objective_matrix_coefficent_values,
+  const cuopt_int_t* constraint_matrix_row_offsets,
+  const cuopt_int_t* constraint_matrix_column_indices,
+  const cuopt_float_t* constraint_matrix_coefficent_values,
+  const char* constraint_sense,
+  const cuopt_float_t* rhs,
+  const cuopt_float_t* lower_bounds,
+  const cuopt_float_t* upper_bounds,
+  cuOptOptimizationProblem* problem_ptr);
+
+/** @brief Create an optimization problem of the form *
+ * @verbatim
+ *                minimize/maximize  c^T x + x^T Q x + offset
+ *                  subject to       bl <= A*x <= bu
+ *                                   l <= x <= u
+ * @endverbatim
+ *
+ * @param[in] num_constraints - The number of constraints.
+ *
+ * @param[in] num_variables - The number of variables.
+ *
+ * @param[in] objective_sense - The objective sense (CUOPT_MINIMIZE for
+ *  minimization or CUOPT_MAXIMIZE for maximization)
+ *
+ * @param[in] objective_offset - An offset to add to the linear objective.
+ *
+ * @param[in] objective_coefficients - A pointer to an array of type cuopt_float_t
+ *  of size num_variables containing the coefficients of the linear objective.
+ *
+ * @param[in] quadratic_objective_matrix_row_offsets - A pointer to an array of type
+ *  cuopt_int_t of size num_variables + 1. quadratic_objective_matrix_row_offsets[i] is the
+ *  index of the first non-zero element of the i-th row of the quadratic objective matrix in
+ *  quadratic_objective_matrix_column_indices and quadratic_objective_matrix_coefficent_values.
+ *  This is part of the compressed sparse row representation of the quadratic objective matrix.
+ *
+ * @param[in] quadratic_objective_matrix_column_indices - A pointer to an array of type
+ *  cuopt_int_t of size quadratic_objective_matrix_row_offsets[num_variables] containing
+ *  the column indices of the non-zero elements of the quadratic objective matrix.
+ *  This is part of the compressed sparse row representation of the quadratic objective matrix.
+ *
+ * @param[in] quadratic_objective_matrix_coefficent_values - A pointer to an array of type
+ *  cuopt_float_t of size quadratic_objective_matrix_row_offsets[num_variables] containing
+ *  the values of the non-zero elements of the quadratic objective matrix.
+ *
+ * @param[in] constraint_matrix_row_offsets - A pointer to an array of type
+ *  cuopt_int_t of size num_constraints + 1. constraint_matrix_row_offsets[i] is the
+ *  index of the first non-zero element of the i-th constraint in
+ *  constraint_matrix_column_indices and constraint_matrix_coefficients.
+ *
+ * @param[in] constraint_matrix_column_indices - A pointer to an array of type
+ *  cuopt_int_t of size constraint_matrix_row_offsets[num_constraints] containing
+ *  the column indices of the non-zero elements of the constraint matrix.
+ *
+ * @param[in] constraint_matrix_coefficients - A pointer to an array of type
+ *  cuopt_float_t of size constraint_matrix_row_offsets[num_constraints] containing
+ *  the values of the non-zero elements of the constraint matrix.
+ *
+ * @param[in] constraint_lower_bounds - A pointer to an array of type
+ *  cuopt_float_t of size num_constraints containing the lower bounds of the constraints.
+ *
+ * @param[in] constraint_upper_bounds - A pointer to an array of type
+ *  cuopt_float_t of size num_constraints containing the upper bounds of the constraints.
+ *
+ * @param[in] variable_lower_bounds - A pointer to an array of type
+ *  cuopt_float_t of size num_variables containing the lower bounds of the variables.
+ *
+ * @param[in] variable_upper_bounds - A pointer to an array of type
+ *  cuopt_float_t of size num_variables containing the upper bounds of the variables.
+ *
+ * @param[out] problem_ptr - A pointer to a cuOptOptimizationProblem.
+ * On output the problem will be created and initialized with the provided data.
+ *
+ * @return A status code indicating success or failure.
+ */
+cuopt_int_t cuOptCreateQuadraticRangedProblem(
+  cuopt_int_t num_constraints,
+  cuopt_int_t num_variables,
+  cuopt_int_t objective_sense,
+  cuopt_float_t objective_offset,
+  const cuopt_float_t* objective_coefficients,
+  const cuopt_int_t* quadratic_objective_matrix_row_offsets,
+  const cuopt_int_t* quadratic_objective_matrix_column_indices,
+  const cuopt_float_t* quadratic_objective_matrix_coefficent_values,
+  const cuopt_int_t* constraint_matrix_row_offsets,
+  const cuopt_int_t* constraint_matrix_column_indices,
+  const cuopt_float_t* constraint_matrix_coefficients,
+  const cuopt_float_t* constraint_lower_bounds,
+  const cuopt_float_t* constraint_upper_bounds,
+  const cuopt_float_t* variable_lower_bounds,
+  const cuopt_float_t* variable_upper_bounds,
+  cuOptOptimizationProblem* problem_ptr);
+
 /** @brief Destroy an optimization problem
  *
  * @param[in, out] problem_ptr - A pointer to a cuOptOptimizationProblem. On

cpp/include/cuopt/linear_programming/optimization_problem.hpp

@@ -186,6 +186,48 @@ class optimization_problem_t {
    * @param objective_offset Objective offset value.
    */
   void set_objective_offset(f_t objective_offset);
+
+  /**
+   * @brief Set the quadratic objective matrix (Q) in CSR format for QPS files.
+   *
+   * @note This is used for quadratic programming problems where the objective
+   * function contains quadratic terms: x^T * Q * x + c^T * x
+   *
+   * @param[in] Q_values Values of the CSR representation of the quadratic objective matrix
+   * @param size_values Size of the Q_values array
+   * @param[in] Q_indices Indices of the CSR representation of the quadratic objective matrix
+   * @param size_indices Size of the Q_indices array
+   * @param[in] Q_offsets Offsets of the CSR representation of the quadratic objective matrix
+   * @param size_offsets Size of the Q_offsets array
+   * @param validate_positive_semi_definite Whether to validate if the matrix is positive semi
+   * definite
+   */
+  void set_quadratic_objective_matrix(const f_t* Q_values,
+                                      i_t size_values,
+                                      const i_t* Q_indices,
+                                      i_t size_indices,
+                                      const i_t* Q_offsets,
+                                      i_t size_offsets,
+                                      bool validate_positive_semi_definite = false);
+
+  /**
+   * @brief Get the quadratic objective matrix offsets
+   * @return const reference to the Q_offsets vector
+   */
+  const std::vector<i_t>& get_quadratic_objective_offsets() const;
+
+  /**
+   * @brief Get the quadratic objective matrix indices
+   * @return const reference to the Q_indices vector
+   */
+  const std::vector<i_t>& get_quadratic_objective_indices() const;
+
+  /**
+   * @brief Get the quadratic objective matrix values
+   * @return const reference to the Q_values vector
+   */
+  const std::vector<f_t>& get_quadratic_objective_values() const;
+
   /**
    * @brief Set the variables (x) lower bounds.
    * @note Setting before calling the solver is optional, default value for all
@@ -341,6 +383,8 @@ class optimization_problem_t {
   const std::vector<std::string>& get_variable_names() const;
   const std::vector<std::string>& get_row_names() const;
 
+  bool has_quadratic_objective() const;
+
   /**
    * @brief Gets the device-side view (with raw pointers), for ease of access
    *        inside cuda kernels
@@ -382,6 +426,12 @@ class optimization_problem_t {
   f_t objective_scaling_factor_{1};
   /** offset of the objective function */
   f_t objective_offset_{0};
+
+  /** Quadratic objective matrix in CSR format (for (1/2) * x^T * Q * x term) */
+  std::vector<i_t> Q_offsets_;
+  std::vector<i_t> Q_indices_;
+  std::vector<f_t> Q_values_;
+
   /** lower bounds of the variables (primal part) */
   rmm::device_uvector<f_t> variable_lower_bounds_;
   /** upper bounds of the variables (primal part) */

cpp/libmps_parser/include/mps_parser/mps_data_model.hpp

@@ -115,6 +115,7 @@ class mps_data_model_t {
    * @param objective_offset Objective offset value.
    */
   void set_objective_offset(f_t objective_offset);
+
   /**
    * @brief Set the variables (x) lower bounds.
    * @note Setting before calling the solver is optional, default value for all is 0.
@@ -356,7 +357,7 @@ class mps_data_model_t {
 
   // QPS-specific data members for quadratic programming support
   /** Quadratic objective matrix in CSR format (for (1/2) * x^T * Q * x term) */
-  std::vector<f_t> Q_objective_;
+  std::vector<f_t> Q_objective_values_;
   std::vector<i_t> Q_objective_indices_;
   std::vector<i_t> Q_objective_offsets_;
 

cpp/libmps_parser/src/mps_data_model.cpp

@@ -201,8 +201,8 @@ void mps_data_model_t<i_t, f_t>::set_quadratic_objective_matrix(const f_t* Q_val
     mps_parser_expects(
       Q_values != nullptr, error_type_t::ValidationError, "Q_values cannot be null");
   }
-  Q_objective_.resize(size_values);
-  std::copy(Q_values, Q_values + size_values, Q_objective_.data());
+  Q_objective_values_.resize(size_values);
+  std::copy(Q_values, Q_values + size_values, Q_objective_values_.data());
 
   if (size_indices != 0) {
     mps_parser_expects(
@@ -421,13 +421,13 @@ i_t mps_data_model_t<i_t, f_t>::get_nnz() const
 template <typename i_t, typename f_t>
 const std::vector<f_t>& mps_data_model_t<i_t, f_t>::get_quadratic_objective_values() const
 {
-  return Q_objective_;
+  return Q_objective_values_;
 }
 
 template <typename i_t, typename f_t>
 std::vector<f_t>& mps_data_model_t<i_t, f_t>::get_quadratic_objective_values()
 {
-  return Q_objective_;
+  return Q_objective_values_;
 }
 
 template <typename i_t, typename f_t>
@@ -457,7 +457,7 @@ std::vector<i_t>& mps_data_model_t<i_t, f_t>::get_quadratic_objective_offsets()
 template <typename i_t, typename f_t>
 bool mps_data_model_t<i_t, f_t>::has_quadratic_objective() const noexcept
 {
-  return !Q_objective_.empty();
+  return !Q_objective_values_.empty();
 }
 
 // NOTE: Explicitly instantiate all types here in order to avoid linker error