R/254_reductions_solvers_qp_solvers_xpress_qpif.R

In CVXR: Disciplined Convex Optimization

#####
## DO NOT EDIT THIS FILE!! EDIT THE SOURCE INSTEAD: rsrc_tree/reductions/solvers/qp_solvers/xpress_qpif.R
#####

## CVXPY SOURCE: reductions/solvers/qp_solvers/xpress_qpif.py
## XPRESS QP solver interface for QP/MIQP problems
##
## Uses QpSolver.apply() for sign-correct A_eq/b_eq/F_ineq/g_ineq data.
## MIP_CAPABLE = TRUE (supports MIQP via columntypes).
##
## Key differences from Python xpress_qpif.py:
##   - R xpress uses xprs_loadproblemdata() convenience function
##   - Controls set via getcontrolinfo() type dispatch
##   - Qobj = symmetrized P (not upper-tri extraction; R convenience handles it)
##   - Duals via getduals()$duals (negated, matching CVXPY)
##   - Objective: LPOBJVAL for LP/QP, MIPOBJVAL for MIP


## XPRESS_STATUS_MAP and .xpress_set_controls are defined in xpress_conif.R
## (loads first in the conic_solvers/ directory).

# -- XPRESS_QP_Solver class ---------------------------------------------------
## CVXPY SOURCE: xpress_qpif.py

#' @keywords internal
XPRESS_QP_Solver <- new_class("XPRESS_QP_Solver", parent = QpSolver,
  package = "CVXR",
  constructor = function() {
    new_object(S7_object(),
      .cache = new.env(parent = emptyenv()),
      MIP_CAPABLE = TRUE,
      BOUNDED_VARIABLES = FALSE,
      SUPPORTED_CONSTRAINTS = list(Zero, NonNeg),
      REQUIRES_CONSTR = FALSE
    )
  }
)

method(solver_name, XPRESS_QP_Solver) <- function(x) XPRESS_SOLVER


# -- solve_via_data ------------------------------------------------------------
## CVXPY SOURCE: xpress_qpif.py solve_via_data
## Receives QP data from QpSolver.apply(): P, q, A_eq, b_eq, F_ineq, g_ineq

method(solve_via_data, XPRESS_QP_Solver) <- function(x, data, warm_start = FALSE,
                                                      verbose = FALSE,
                                                      solver_opts = list(), ...) {
  if (!requireNamespace("xpress", quietly = TRUE))
    cli_abort("Package {.pkg xpress} is required but not installed.")

  q_vec <- data[["q"]]
  nvars <- length(q_vec)

  A_eq    <- data[["A_eq"]]
  b_eq    <- data[["b_eq"]]
  F_ineq  <- data[["F_ineq"]]
  g_ineq  <- data[["g_ineq"]]

  len_eq   <- nrow(A_eq)
  len_ineq <- nrow(F_ineq)

  ## MIP detection
  bool_idx <- data[["bool_idx"]] %||% integer(0)
  int_idx  <- data[["int_idx"]]  %||% integer(0)
  is_mip   <- length(bool_idx) > 0L || length(int_idx) > 0L

  ## Variable types (R 1-based indices into character vector)
  columntypes <- rep("C", nvars)
  if (length(bool_idx) > 0L) columntypes[bool_idx] <- "B"
  if (length(int_idx) > 0L)  columntypes[int_idx]  <- "I"

  ## Variable bounds
  lb <- rep(-Inf, nvars)
  ub <- rep(Inf, nvars)
  if (length(bool_idx) > 0L) {
    lb[bool_idx] <- pmax(lb[bool_idx], 0)
    ub[bool_idx] <- pmin(ub[bool_idx], 1)
  }

  ## -- Build problem data list for xprs_loadproblemdata -----------------------
  ## Load equality constraints + quadratic objective at once
  problemdata <- list(probname = "CVX_xpress_qp", objcoef = q_vec)

  ## Quadratic objective: symmetrize P, keep sparse
  P <- data[[SD_P]]
  if (!is.null(P) && length(P@x) > 0L) {
    P_sym <- (P + Matrix::t(P)) / 2
    if (!inherits(P_sym, "sparseMatrix"))
      P_sym <- methods::as(P_sym, "CsparseMatrix")
    problemdata$Qobj <- P_sym
  }

  ## Equality constraints (always pass A — xprs_loadproblemdata requires it)
  if (len_eq > 0L) {
    problemdata$A <- A_eq
    problemdata$rhs <- b_eq
    problemdata$rowtype <- rep("E", len_eq)
  } else {
    problemdata$A <- Matrix::sparseMatrix(i = integer(0), j = integer(0),
                                           x = numeric(0), dims = c(0L, nvars))
    problemdata$rhs <- numeric(0)
    problemdata$rowtype <- character(0)
  }

  problemdata$lb <- lb
  problemdata$ub <- ub
  if (is_mip) problemdata$columntypes <- columntypes

  ## -- Create problem and load ------------------------------------------------
  prob <- xpress::createprob()
  on.exit(xpress::destroyprob(prob), add = TRUE)

  xpress::xprs_loadproblemdata(prob, problemdata = problemdata)

  ## -- Add inequality constraints ---------------------------------------------
  if (len_ineq > 0L) {
    F_csr <- methods::as(F_ineq, "RsparseMatrix")
    xpress::addrows(prob,
      rowtype = rep("L", len_ineq),
      rhs     = g_ineq,
      start   = F_csr@p,
      colind  = F_csr@j,
      rowcoef = F_csr@x
    )
  }

  ## -- Set solver controls ----------------------------------------------------
  ## Defaults matching CVXPY
  defaults <- list()
  if (!"BARGAPTARGET" %in% toupper(names(solver_opts)))
    defaults[["BARGAPTARGET"]] <- 1e-30
  if (!"FEASTOL" %in% toupper(names(solver_opts)))
    defaults[["FEASTOL"]] <- 1e-9
  .xpress_set_controls(prob, defaults)

  ## Verbosity
  if (verbose) {
    .xpress_set_controls(prob, list(OUTPUTLOG = 1L, MIPLOG = 2L, LPLOG = 1L))
  } else {
    .xpress_set_controls(prob, list(OUTPUTLOG = 0L, MIPLOG = 0L, LPLOG = 0L))
  }

  ## User-specified controls (skip write_mps)
  user_opts <- solver_opts[!names(solver_opts) %in% c("write_mps")]
  if (length(user_opts) > 0L) .xpress_set_controls(prob, user_opts)

  ## Write MPS if requested
  if (!is.null(solver_opts[["write_mps"]])) {
    tryCatch(xpress::writeprob(prob, solver_opts[["write_mps"]]),
             error = function(e) NULL)
  }

  ## -- MIP warm-start ---------------------------------------------------------
  if (warm_start && is_mip) {
    mip_indices <- c(bool_idx, int_idx) - 1L  ## convert to 0-based
    dots <- list(...)
    solver_cache <- dots[["solver_cache"]]
    if (!is.null(solver_cache) && exists(XPRESS_SOLVER, envir = solver_cache)) {
      cached <- get(XPRESS_SOLVER, envir = solver_cache)
      if (!is.null(cached$x) && length(cached$x) == nvars) {
        warm_vals <- cached$x[c(bool_idx, int_idx)]
        tryCatch(
          xpress::addmipsol(prob, warm_vals, mip_indices, "warmstart"),
          error = function(e) NULL
        )
      }
    }
  }

  ## -- Optimize ---------------------------------------------------------------
  tryCatch(
    xpress::xprs_optimize(prob),
    error = function(e) NULL
  )

  ## -- Extract results --------------------------------------------------------
  result <- list()

  sol_status <- xpress::getintattrib(prob, xpress:::SOLSTATUS)
  result$status <- as.character(sol_status)
  result$.len_eq <- len_eq
  result$.is_mip <- is_mip

  status <- XPRESS_STATUS_MAP[[result$status]]
  if (is.null(status)) status <- SOLVER_ERROR

  if (status %in% SOLUTION_PRESENT) {
    result$x <- xpress::getsolution(prob)$x
    result$objval <- if (is_mip) {
      xpress::getdblattrib(prob, xpress:::MIPOBJVAL)
    } else {
      xpress::getdblattrib(prob, xpress:::LPOBJVAL)
    }

    ## Duals (LP/QP only, not MIP)
    nrows <- len_eq + len_ineq
    if (!is_mip && nrows > 0L) {
      result$duals <- xpress::getduals(prob, first = 0L,
                                        last = nrows - 1L)$duals
    }
  }

  tryCatch(
    result$solve_time <- xpress::getdblattrib(prob, xpress:::TIME),
    error = function(e) NULL
  )

  ## Cache for future warm-starts
  dots <- list(...)
  solver_cache <- dots[["solver_cache"]]
  if (!is.null(solver_cache)) {
    assign(XPRESS_SOLVER, result, envir = solver_cache)
  }

  result
}

# -- reduction_invert ----------------------------------------------------------
## CVXPY SOURCE: xpress_qpif.py invert
## Dual sign: negate ALL duals (matching CVXPY xpress_qpif.py)

method(reduction_invert, XPRESS_QP_Solver) <- function(x, solution, inverse_data, ...) {
  attr_list <- list()

  status <- XPRESS_STATUS_MAP[[solution$status]]
  if (is.null(status)) status <- SOLVER_ERROR

  ## Timing
  if (!is.null(solution$solve_time))
    attr_list[[RK_SOLVE_TIME]] <- solution$solve_time

  if (status %in% SOLUTION_PRESENT) {
    opt_val <- solution$objval + inverse_data[[SD_OFFSET]]

    ## Primal variables
    primal_vars <- list()
    primal_vars[[as.character(inverse_data[[SOLVER_VAR_ID]])]] <- solution$x

    ## Dual variables: negate ALL duals (CVXPY: y = -np.array(getDual()))
    is_mip <- isTRUE(solution$.is_mip)

    if (!is_mip && !is.null(solution$duals)) {
      y <- -solution$duals
      len_eq <- solution$.len_eq

      eq_dual <- if (len_eq > 0L) {
        get_dual_values(
          y[seq_len(len_eq)],
          extract_dual_value,
          inverse_data[[SOLVER_EQ_CONSTR]]
        )
      } else {
        list()
      }

      ineq_dual <- if (len_eq < length(y)) {
        get_dual_values(
          y[(len_eq + 1L):length(y)],
          extract_dual_value,
          inverse_data[[SOLVER_NEQ_CONSTR]]
        )
      } else {
        list()
      }

      dual_vars <- c(eq_dual, ineq_dual)
    } else {
      dual_vars <- list()
    }

    Solution(status, opt_val, primal_vars, dual_vars, attr_list)
  } else {
    failure_solution(status, attr_list)
  }
}

# -- print ---------------------------------------------------------------------

method(print, XPRESS_QP_Solver) <- function(x, ...) {
  cat("XPRESS_QP_Solver()\n")
  invisible(x)
}