# Integer Programming

Consider the following optimization problem.

$\begin{array}{ll} \mbox{Maximize} & x_1 + 2x_2 - 0.1x_3 - 3x_4\\ \mbox{subject to} & x_1, x_2, x_3, x_4 >= 0\\ & x_1 + x_2 <= 5\\ & 2x_1 - x_2 >= 0\\ & -x_1 + 3x_2 >= 0\\ & x_3 + x_4 >= 0.5\\ & x_3 >= 1.1\\ & x_3 \mbox{ is integer.} \end{array}$

CVXR provides constructors for the integer and boolean variables via the parameter integer = TRUE or boolean = TRUE to the Variable() function. These can be combined with vstack (analog of rbind) or hstack (analog of cbind) to construct new expressions.

The above problem now in CVXR.

y1 <- Variable(2)
y2 <- Variable(1, integer = TRUE)
y3 <- Variable(1)
x <- vstack(y1, y2, y3) ## Create x expression
C <- matrix(c(1, 2, -0.1, -3), nrow = 1)
objective <- Maximize(C %*% x)
constraints <- list(
x >= 0,
x[1] + x[2] <= 5,
2 * x[1] - x[2] >= 0,
-x[1] + 3 * x[2] >= 0,
x[3] + x[4] >= 0.5,
x[3] >= 1.1)
problem <- Problem(objective, constraints)

We can solve this problem as usual using the default ECOS (actually ECOS_BB) solver and obtain the optimal value as well as the solution.

result <- solve(problem, solver = "ECOS_BB")
cat(sprintf("Optimal value: %.3f\n", result$value)) ## Optimal value: 8.133 ecos_solution <- result$getValue(x)

## Alternative Solvers

We can try other solvers and compare the solutions obtained, like LPSOLVE and GLPK provided the respective R packages are installed as documented in the tutorial Using Other Solvers.

Note: LPSOLVE is now orphaned on CRAN and so no longer supported.

##result <- solve(problem, solver = "LPSOLVE")
##lpsolve_solution <- result$getValue(x) result <- solve(problem, solver = "GLPK") glpk_solution <- result$getValue(x)

We can also try a commercial solver, Gurobi, that can handle integer programs. This requires installation of the Gurobi solver, discussed in Using Other Solvers.

result <- solve(problem, solver = "GUROBI")
gurobi_solution <- result$getValue(x) Finally, two other commercial solvers, MOSEK and CPLEX. result <- solve(problem, solver = "CPLEX") cplex_solution <- result$getValue(x)
result <- solve(problem, solver = "MOSEK")
mosek_solution <- result$getValue(x) Below is the table is solutions from all the solvers we used. solutions <- data.frame(ECOS = ecos_solution, ## LPSOLVE = lpsolve_solution, GLPK = glpk_solution, GUROBI = gurobi_solution, MOSEK = mosek_solution, CPLEX = cplex_solution) row.names(solutions) <- c("$x_1$", "$x_2$", "$x_3$", "$x_4$") knitr::kable(solutions, format = "html") %>% kable_styling("striped") %>% column_spec(1:4, background = "#ececec") ECOS GLPK GUROBI MOSEK CPLEX $$x_1$$ 1.6666703 1.666667 1.666667 1.666667 1.666667 $$x_2$$ 3.3333291 3.333333 3.333333 3.333333 3.333333 $$x_3$$ 2.0000704 2.000000 2.000000 2.000000 2.000000 $$x_4$$ -0.0000018 0.000000 0.000000 0.000000 0.000000 ## Office Assignment Problem For a slightly more involved example, we consider the office assignment problem. The goal is to assign six people, Marcelo, Rakesh, Peter, Tom, Marjorie, and Mary Ann, to seven offices. Each office can have no more than one person, and each person gets exactly one office. So there will be one empty office. People can give preferences for the offices, and their preferences are considered based on their seniority. Some offices have windows, some do not, and one window is smaller than others. Additionally, Peter and Tom often work together, so should be in adjacent offices. Marcelo and Rakesh often work together, and should be in adjacent offices. draw_office_layout() The office layout is shown above. Offices 1, 2, 3, and 4 are inside offices (no windows). Offices 5, 6, and 7 have windows, but the window in office 5 is smaller than the other two. We begin by recording the names of the people and offices. people <- c('Mary Ann', 'Marjorie', 'Tom', 'Peter', 'Marcelo', 'Rakesh') offices <- c('Office 1', 'Office 2', 'Office 3', 'Office 4','Office 5', 'Office 6', 'Office 7') We also have the office preferences of each person for each of the seven offices along with seniority data which is used to scale the office preferences. preference_matrix <- matrix( c(0, 0, 0, 0, 10, 40, 50, 0, 0, 0, 0, 20, 40, 40, 0, 0, 0, 0, 30, 40, 30, 1, 3, 3, 3, 10, 40, 40, 3, 4, 1, 2, 10, 40, 40, 10, 10, 10, 10, 20, 20, 20), byrow = TRUE, nrow = length(people)) rownames(preference_matrix) <- people colnames(preference_matrix) <- offices seniority <- c(9, 10, 5, 3, 1.5, 2) weightvector <- seniority / sum(seniority) PM <- diag(weightvector) %*% preference_matrix We define the the occupancy variable which indicates, using values 1 or 0, who occupies which office. occupy <- Variable(length(people), length(offices), integer = TRUE) The objective is to maximize the satisfaction of the preferences weighted by seniority constrained by the fact the a person can only occupy a single office and no office can have more than 1 person. objective <- Maximize(sum_entries(multiply(PM, occupy))) constraints <- list( occupy >= 0, occupy <= 1, sum_entries(occupy, axis = 1) == 1, sum_entries(occupy, axis = 2) <= 1 ) We further add the constraint that Tom (person 3) and Peter (person 4) should be no more than one office away, and ditto for Marcelo (person 5) and Rakesh (person 6). tom_peter <- list( occupy[3, 1] + sum_entries(occupy[4, ]) - occupy[4, 2] <= 1, occupy[3, 2] + sum_entries(occupy[4, ]) - occupy[4, 1] - occupy[4, 3] - occupy[4, 5] <= 1, occupy[3, 3] + sum_entries(occupy[4, ]) - occupy[4, 2] - occupy[4, 4] - occupy[4, 6] <= 1, occupy[3, 4] + sum_entries(occupy[4, ]) - occupy[4, 3] - occupy[4, 7] <= 1, occupy[3, 5] + sum_entries(occupy[4, ]) - occupy[4, 2] - occupy[4, 6] <= 1, occupy[3, 6] + sum_entries(occupy[4, ]) - occupy[4, 3] - occupy[4, 5] - occupy[4, 7] <= 1, occupy[3, 7] + sum_entries(occupy[4, ]) - occupy[4, 4] - occupy[4, 6] <= 1 ) marcelo_rakesh <- list( occupy[5, 1] + sum_entries(occupy[6, ]) - occupy[6, 2] <= 1, occupy[5, 2] + sum_entries(occupy[6, ]) - occupy[6, 1] - occupy[6, 3] - occupy[6, 5] <= 1, occupy[5, 3] + sum_entries(occupy[6, ]) - occupy[6, 2] - occupy[6, 4] - occupy[6, 6] <= 1, occupy[5, 4] + sum_entries(occupy[6, ]) - occupy[6, 3] - occupy[6, 7] <= 1, occupy[5, 5] + sum_entries(occupy[6, ]) - occupy[6, 2] - occupy[6, 6] <= 1, occupy[5, 6] + sum_entries(occupy[6, ]) - occupy[6, 3] - occupy[6, 5] - occupy[6, 7] <= 1, occupy[5, 7] + sum_entries(occupy[6, ]) - occupy[6, 4] - occupy[6, 6] <= 1 ) constraints <- c(constraints, tom_peter, marcelo_rakesh) We are now ready to solve the problem. problem <- Problem(objective, constraints) ecos_result <- solve(problem, solver = "ECOS_BB") ecos_soln <- round(ecos_result$getValue(occupy), 0)
rownames(ecos_soln) <- people
colnames(ecos_soln) <- offices

We are now ready to plot the solution (after accounting for the WC).

office_assignment <- apply(ecos_soln, 1, which.max)
office_occupants <- names(office_assignment)[match(c(5:7, 1:4), office_assignment)]
office_occupants[is.na(office_occupants)] <- "Empty"
draw_office_layout(c("WC", office_occupants))

## Session Info

sessionInfo()
## R version 4.2.1 (2022-06-23)
## Platform: x86_64-apple-darwin21.6.0 (64-bit)
## Running under: macOS Ventura 13.0
##
## Matrix products: default
## BLAS:   /usr/local/Cellar/openblas/0.3.21/lib/libopenblasp-r0.3.21.dylib
## LAPACK: /usr/local/Cellar/r/4.2.1_4/lib/R/lib/libRlapack.dylib
##
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
##
## attached base packages:
## [1] stats     graphics  grDevices datasets  utils     methods   base
##
## other attached packages:
## [1] ggplot2_3.3.6    kableExtra_1.3.4 CVXR_1.0-11
##
## loaded via a namespace (and not attached):
##  [1] Rcpp_1.0.9        svglite_2.1.0     gurobi_9.5-2      lattice_0.20-45
##  [5] rcbc_0.1.0.9001   Rmosek_10.0.25    assertthat_0.2.1  digest_0.6.30
##  [9] utf8_1.2.2        gmp_0.6-6         slam_0.1-50       R6_2.5.1
## [13] evaluate_0.17     httr_1.4.4        highr_0.9         blogdown_1.13
## [17] pillar_1.8.1      rlang_1.0.6       rstudioapi_0.14   jquerylib_0.1.4
## [21] Matrix_1.5-1      rmarkdown_2.17    labeling_0.4.2    webshot_0.5.4
## [25] stringr_1.4.1     bit_4.0.4         munsell_0.5.0     compiler_4.2.1
## [29] xfun_0.34         pkgconfig_2.0.3   systemfonts_1.0.4 Rglpk_0.6-4
## [33] htmltools_0.5.3   Rcplex_0.3-5      tidyselect_1.2.0  tibble_3.1.8
## [37] bookdown_0.29     codetools_0.2-18  fansi_1.0.3       viridisLite_0.4.1
## [41] dplyr_1.0.10      withr_2.5.0       grid_4.2.1        jsonlite_1.8.3
## [45] gtable_0.3.1      lifecycle_1.0.3   DBI_1.1.3         magrittr_2.0.3
## [49] scales_1.2.1      cli_3.4.1         stringi_1.7.8     cachem_1.0.6
## [53] farver_2.1.1      Rmpfr_0.8-9       xml2_1.3.3        bslib_0.4.0
## [57] generics_0.1.3    vctrs_0.5.0       tools_4.2.1       bit64_4.0.5
## [61] glue_1.6.2        fastmap_1.1.0     yaml_2.3.6        colorspace_2.0-3
## [65] cccp_0.2-9        rvest_1.0.3       ECOSolveR_0.5.4   knitr_1.40
## [69] sass_0.4.2

R Markdown