Skip to content

Write quadratic terms to mps file #949

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
rg20 wants to merge 11 commits into
base: release/26.04
Choose a base branch
from
Open

Write quadratic terms to mps file #949

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
11 commits
Select commit Hold shift + click to select a range
85902c3
speedup symmetrization
rg20 Feb 10, 2026
bbdbbdb
separate out heper function
rg20 Feb 10, 2026
e9f7a79
Support quadobj writing to mps
rg20 Feb 12, 2026
5b5a2bb
add missing header
rg20 Feb 13, 2026
eb79da4
Use kahan summation for post solve residual computation
rg20 Mar 5, 2026
1b94b55
Handle core review comment
rg20 Mar 16, 2026
0f3ca8e
Merge remote-tracking branch 'upstream/release/26.04' into qp_improve…
rg20 Mar 20, 2026
2c05667
Remove debug print
rg20 Mar 23, 2026
172ea69
Merge remote-tracking branch 'upstream/release/26.04' into qp_improve…
rg20 Mar 23, 2026
55aa43a
Remove redundant function
rg20 Mar 23, 2026
f6a8240
assert guard for workspace
rg20 Mar 23, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
1 change: 1 addition & 0 deletions cpp/libmps_parser/CMakeLists.txt
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -111,6 +111,7 @@ add_library(cuopt::mps_parser ALIAS mps_parser)
target_include_directories(mps_parser
PRIVATE
"${CMAKE_CURRENT_SOURCE_DIR}/../thirdparty"
"${CMAKE_CURRENT_SOURCE_DIR}/../src"
"${CMAKE_CURRENT_SOURCE_DIR}/src"
PUBLIC
"$<BUILD_INTERFACE:${CMAKE_CURRENT_SOURCE_DIR}/include>"
Expand Down
18 changes: 16 additions & 2 deletions cpp/libmps_parser/include/mps_parser/mps_writer.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,16 +1,18 @@
/* clang-format off */
/*
* SPDX-FileCopyrightText: Copyright (c) 2025, NVIDIA CORPORATION & AFFILIATES. All rights reserved.
* SPDX-FileCopyrightText: Copyright (c) 2025-2026, NVIDIA CORPORATION & AFFILIATES. All rights reserved.
* SPDX-License-Identifier: Apache-2.0
*/
/* clang-format on */

#pragma once

#include <mps_parser/data_model_view.hpp>
#include <mps_parser/mps_data_model.hpp>

#include <stdarg.h>
#include <limits>
#include <memory>
#include <string>
#include <unordered_map>
#include <unordered_set>
Expand All @@ -31,10 +33,16 @@ class mps_writer_t {
* @brief Ctor. Takes a data model view as input and writes it out as a MPS formatted file
*
* @param[in] problem Data model view to write
* @param[in] file Path to the MPS file to write
*/
mps_writer_t(const data_model_view_t<i_t, f_t>& problem);

/**
* @brief Ctor. Takes a data model as input and writes it out as a MPS formatted file
*
* @param[in] problem Data model to write
*/
mps_writer_t(const mps_data_model_t<i_t, f_t>& problem);

/**
* @brief Writes the problem to an MPS formatted file
*
Expand All @@ -43,7 +51,13 @@ class mps_writer_t {
void write(const std::string& mps_file_path);

private:
// Owned view (created when constructing from mps_data_model_t)
std::unique_ptr<data_model_view_t<i_t, f_t>> owned_view_;
// Reference to the view (either external or owned)
const data_model_view_t<i_t, f_t>& problem_;

// Helper to create view from data model
static data_model_view_t<i_t, f_t> create_view(const mps_data_model_t<i_t, f_t>& model);
}; // class mps_writer_t

} // namespace cuopt::mps_parser
140 changes: 140 additions & 0 deletions cpp/libmps_parser/src/mps_writer.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -8,14 +8,17 @@
#include <mps_parser/mps_writer.hpp>

#include <mps_parser/data_model_view.hpp>
#include <mps_parser/mps_data_model.hpp>
#include <utilities/error.hpp>
#include <utilities/sparse_matrix_helpers.hpp>

#include <cmath>
#include <fstream>
#include <iomanip>
#include <iostream>
#include <limits>
#include <map>
#include <memory>

namespace cuopt::mps_parser {

Expand All @@ -24,6 +27,92 @@ mps_writer_t<i_t, f_t>::mps_writer_t(const data_model_view_t<i_t, f_t>& problem)
{
}

template <typename i_t, typename f_t>
data_model_view_t<i_t, f_t> mps_writer_t<i_t, f_t>::create_view(
const mps_data_model_t<i_t, f_t>& model)
{
data_model_view_t<i_t, f_t> view;

// Set basic data
view.set_maximize(model.get_sense());

// Constraint matrix
const auto& A_values = model.get_constraint_matrix_values();
const auto& A_indices = model.get_constraint_matrix_indices();
const auto& A_offsets = model.get_constraint_matrix_offsets();
if (!A_values.empty()) {
view.set_csr_constraint_matrix(A_values.data(),
static_cast<i_t>(A_values.size()),
A_indices.data(),
static_cast<i_t>(A_indices.size()),
A_offsets.data(),
static_cast<i_t>(A_offsets.size()));
}
Comment on lines +43 to +50
Copy link
Copy Markdown

@coderabbitai coderabbitai bot Mar 16, 2026

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

⚠️ Potential issue | 🔴 Critical

Don’t skip CSR setup when A_values is empty

Line 43 should not gate on values. A valid zero-nnz CSR still needs offsets; otherwise write() can read constraint_matrix_offsets[row_id + 1] out of bounds for constrained-but-empty matrices.

✅ Proposed fix 
-  if (!A_values.empty()) {
+  if (!A_offsets.empty()) {
     view.set_csr_constraint_matrix(A_values.data(),
                                    static_cast<i_t>(A_values.size()),
                                    A_indices.data(),
                                    static_cast<i_t>(A_indices.size()),
                                    A_offsets.data(),
                                    static_cast<i_t>(A_offsets.size()));
   }
As per coding guidelines, "Validate algorithm correctness in optimization logic: simplex pivots, branch-and-bound decisions, routing heuristics, and constraint/objective handling must produce correct results."
🤖 Prompt for AI Agents 
Verify each finding against the current code and only fix it if needed.

In `@cpp/libmps_parser/src/mps_writer.cpp` around lines 43 - 50, The code
currently skips calling view.set_csr_constraint_matrix when A_values.empty(),
which breaks zero-nnz CSR matrices because write() later accesses
constraint_matrix_offsets[row_id + 1]; instead ensure set_csr_constraint_matrix
is invoked even when A_values is empty by removing the A_values.empty() guard
and calling view.set_csr_constraint_matrix with the existing A_offsets (and
A_values/A_indices pointers with size 0 as appropriate) whenever the constraint
matrix structure (A_offsets) is present or the matrix is constrained; update the
call site (the block around view.set_csr_constraint_matrix) to pass correct zero
sizes (static_cast<i_t>(0)) for values/indices when A_values/A_indices are empty
so constraint_matrix_offsets is always initialized for write().


// Constraint bounds
const auto& b = model.get_constraint_bounds();
if (!b.empty()) { view.set_constraint_bounds(b.data(), static_cast<i_t>(b.size())); }

// Objective coefficients
const auto& c = model.get_objective_coefficients();
if (!c.empty()) { view.set_objective_coefficients(c.data(), static_cast<i_t>(c.size())); }

view.set_objective_scaling_factor(model.get_objective_scaling_factor());
view.set_objective_offset(model.get_objective_offset());

// Variable bounds
const auto& lb = model.get_variable_lower_bounds();
const auto& ub = model.get_variable_upper_bounds();
if (!lb.empty()) { view.set_variable_lower_bounds(lb.data(), static_cast<i_t>(lb.size())); }
if (!ub.empty()) { view.set_variable_upper_bounds(ub.data(), static_cast<i_t>(ub.size())); }

// Variable types
const auto& var_types = model.get_variable_types();
if (!var_types.empty()) {
view.set_variable_types(var_types.data(), static_cast<i_t>(var_types.size()));
}

// Row types
const auto& row_types = model.get_row_types();
if (!row_types.empty()) {
view.set_row_types(row_types.data(), static_cast<i_t>(row_types.size()));
}

// Constraint bounds (lower/upper)
const auto& cl = model.get_constraint_lower_bounds();
const auto& cu = model.get_constraint_upper_bounds();
if (!cl.empty()) { view.set_constraint_lower_bounds(cl.data(), static_cast<i_t>(cl.size())); }
if (!cu.empty()) { view.set_constraint_upper_bounds(cu.data(), static_cast<i_t>(cu.size())); }

// Names
view.set_problem_name(model.get_problem_name());
view.set_objective_name(model.get_objective_name());
view.set_variable_names(model.get_variable_names());
view.set_row_names(model.get_row_names());

// Quadratic objective
const auto& Q_values = model.get_quadratic_objective_values();
const auto& Q_indices = model.get_quadratic_objective_indices();
const auto& Q_offsets = model.get_quadratic_objective_offsets();
if (!Q_values.empty()) {
view.set_quadratic_objective_matrix(Q_values.data(),
static_cast<i_t>(Q_values.size()),
Q_indices.data(),
static_cast<i_t>(Q_indices.size()),
Q_offsets.data(),
static_cast<i_t>(Q_offsets.size()));
}

return view;
}

template <typename i_t, typename f_t>
mps_writer_t<i_t, f_t>::mps_writer_t(const mps_data_model_t<i_t, f_t>& problem)
: owned_view_(std::make_unique<data_model_view_t<i_t, f_t>>(create_view(problem))),
problem_(*owned_view_)
{
}

template <typename i_t, typename f_t>
void mps_writer_t<i_t, f_t>::write(const std::string& mps_file_path)
{
Expand Down Expand Up @@ -280,6 +369,57 @@ void mps_writer_t<i_t, f_t>::write(const std::string& mps_file_path)
}
}

// QUADOBJ section for quadratic objective terms (if present)
// MPS format: QUADOBJ stores upper triangular elements (row <= col)
// MPS uses (1/2) x^T H x, cuOpt uses x^T Q x
// For equivalence: H[i,j] = Q[i,j] + Q[j,i] (works for both diagonal and off-diagonal)
// We symmetrize Q first (H = Q + Q^T), then extract upper triangular
if (problem_.has_quadratic_objective()) {
auto Q_values_span = problem_.get_quadratic_objective_values();
auto Q_indices_span = problem_.get_quadratic_objective_indices();
auto Q_offsets_span = problem_.get_quadratic_objective_offsets();

// Copy span data to local vectors for indexed access
std::vector<f_t> Q_values(Q_values_span.data(), Q_values_span.data() + Q_values_span.size());
std::vector<i_t> Q_indices(Q_indices_span.data(),
Q_indices_span.data() + Q_indices_span.size());
std::vector<i_t> Q_offsets(Q_offsets_span.data(),
Q_offsets_span.data() + Q_offsets_span.size());

if (Q_values.size() > 0) {
// Symmetrize Q: compute H = Q + Q^T
std::vector<f_t> H_values;
std::vector<i_t> H_indices;
std::vector<i_t> H_offsets;
cuopt::symmetrize_csr<i_t, f_t>(
Q_values, Q_indices, Q_offsets, H_values, H_indices, H_offsets);

i_t n_rows = static_cast<i_t>(H_offsets.size()) - 1;

mps_file << "QUADOBJ\n";

// Write upper triangular entries from symmetric H
for (i_t i = 0; i < n_rows; ++i) {
std::string row_name = static_cast<size_t>(i) < problem_.get_variable_names().size()
? problem_.get_variable_names()[i]
: "C" + std::to_string(i);

for (i_t p = H_offsets[i]; p < H_offsets[i + 1]; ++p) {
i_t j = H_indices[p];
f_t v = H_values[p];

// Only write upper triangular (i <= j)
if (i <= j && v != f_t(0)) {
std::string col_name = static_cast<size_t>(j) < problem_.get_variable_names().size()
? problem_.get_variable_names()[j]
: "C" + std::to_string(j);
mps_file << " " << row_name << " " << col_name << " " << v << "\n";
}
}
}
}
}

mps_file << "ENDATA\n";
mps_file.close();
}
Expand Down
Loading
Loading