lpSolve in R with Character and Column Sum Contraints

Question

I have a list of rooms, the rooms' maximum square feet, programs, the programs' maximum square feet, and the values of how well a room matches (match #) the intended program use. With help, I have been able to maximize match # and square feet use for one program per room. However, I would like to take this analysis one step further and allow for multiple programs in the same room or multiples of the same program if it has the highest match #, so long as the multiples still fit within the square foot requirements. Moreover, I would like to tell lpSolve overall that I only want "x" number of Offices, "y" number of Studios, etc. throughout the entire building. Here is my data and code thus far:

program.size <- c(120,320,300,800,500,1000,500,1000,1500,400,1500,2000)
room.size <- c(1414,682,1484,2938,1985,1493,427,1958,708,581,1485,652,727,2556,1634,187,2174,205,1070,2165,1680,1449,1441,2289,986,298,590,2925)
(obj.vals <- matrix(c(3,4,2,8,3,7,4,8,6,4,7,7,
                  3,4,2,8,3,7,4,8,6,4,7,7,
                  4,5,3,7,4,6,5,7,5,3,6,6,
                  2,3,1,7,2,6,3,7,7,5,6,6,
                  4,5,3,7,4,6,5,7,5,3,6,6,
                  3,6,4,8,5,7,4,8,7,7,7,7,
                  3,4,2,8,3,7,4,8,6,4,7,7,
                  4,5,3,7,4,6,5,7,5,3,6,6,
                  6,7,5,7,6,6,7,7,5,3,6,6,
                  6,7,5,7,6,6,7,7,5,3,6,6,
                  5,6,6,6,5,7,8,6,4,2,5,5,
                  6,7,5,7,6,6,7,7,5,3,6,6,
                  6,7,5,7,6,6,7,7,5,3,6,6,
                  3,4,4,8,3,9,6,8,6,4,7,7,
                  3,4,2,6,3,5,4,6,6,4,5,5,
                  4,5,3,5,4,4,5,5,5,3,4,4,
                  5,6,4,8,5,7,6,8,6,4,7,7,
                  5,6,4,8,5,7,6,8,6,4,7,7,
                  4,5,5,7,4,8,7,7,5,3,6,6,
                  5,6,4,8,5,7,6,8,6,4,7,7,
                  3,4,2,6,3,5,4,6,6,4,5,5,
                  5,6,4,8,5,7,6,8,6,4,7,7,
                  5,6,4,8,5,7,6,8,6,4,7,7,
                  5,4,4,6,5,5,6,6,6,6,7,5,
                  6,5,5,5,6,4,5,5,5,7,6,4,
                  4,5,3,7,4,6,5,7,7,5,6,6,
                  6,5,5,5,6,4,5,5,5,7,6,4,
                  3,4,4,6,3,7,6,6,6,4,5,5), nrow=12))
rownames(obj.vals) <- c("Enclosed Offices", "Open Office", "Reception / Greeter", "Studio / Classroom",
                    "Conference / Meeting Room", "Gallery", "Public / Lobby / Waiting", 
                    "Collaborative Space", "Mechanical / Support", "Storage / Archives", 
                    "Fabrication", "Performance")
(obj.adj <- obj.vals * outer(program.size, room.size, "<="))


nr <- nrow(obj.adj)
nc <- ncol(obj.adj)

library(lpSolve)
obj <- as.vector(obj.adj)
con <- t(1*sapply(1:nc, function(x) rep(1:nc == x, each=nr)))
dir <- rep("<=", nc)
rhs <- rep(1, nc)
mod <- lp("max", obj, con, dir, rhs, all.bin=TRUE)

final <- matrix(mod$solution, nrow=nr)

And so now my question is how can I allow the solver to maximize square foot use and match # within each room (column) and allow either multiple of the same programs or a combination of programs to accomplish this? I know I would have to lift the "<= 1" restriction in "mod", but I can't figure out how to let it find the best fit in each room and then ultimately, overall.

The solution that should come for room [,1] is:

$optimum
33

And it's going to try to fit 11 Enclosed Offices within the room which scores a much higher optimum match # than 1 Collaborative Space (8 matches) and 1 Storage / Archives (4 matches) for a total of 12 matches.

And so this leads to my next question about limiting the overall number of certain programs within my solution matrix. I assume it would include some kind of

as.numeric(data$EnclosedOffices "<=" 5)

but I also can't figure out how to limit that. These numbers would be different for all the programs.

Thanks for any and all help and feel free to ask for any clarification.

---

Update: Constraints

• Maximize the match # (obj.vals) for each room.
• Maximize the program.size square feet within each room.size square feet. Do this by either using the same program multiple times (5 x Enclosed Offices) or combining programs (1 Collaborative Space and 1 Performance).
• Restrict and/or force the number of programs returned in the solution. The programs can be split up among rooms so long as it doesn't go over the maximum number I provide for that program in all the rooms. (Only 5 Enclosed Offices, 8 Studios/Classrooms, 1 Fabrication, etc. can be selected across all 28 columns.

Answer 1

If you use the R package lpSolveAPI (a wrapper for lpSolve) then it becomes a little easier. First, take a look at the mathematical formulation (an Integer Program) and then I show you the code to solve your problem.

Formulation

Let X_r_p be the decision variable that takes on positive integer values.

X_r_p = Number of programs of type p assigned to room r (In all your problem will have 28*12=336 decision variables)

Objective Function

Maximize matching score

Max sum(r) sum(p) C_r_p * X_r_p # Here C_r_p is the score of assigning p to room r

Subject to

Room Area Limitation Constraint

Sum(p) Max_area_p * X_r_p <= Room Size (r) for each room r

(We will have 28 such constraints)

Restrict the Number of programs Constraint

Sum(r) X_r_p <= Max_allowable(p) for each program p

(We will have 12 such constraints)

 X_r_p >= 0, Integer

That is the formulation in all. 336 columns and 40 rows.

Implementatin in R

Here's an implementation in R, using lpSolveAPI. Note: Since the OP didn't provide the max_allowable programs in the building, I generated my own data for max_programs.

program.size <- c(120,320,300,800,500,1000,500,1000,1500,400,1500,2000)
    room.size <- c(1414,682,1484,2938,1985,1493,427,1958,708,581,1485,652,727,2556,1634,187,2174,205,1070,2165,1680,1449,1441,2289,986,298,590,2925)
    (obj.vals <- matrix(c(3,4,2,8,3,7,4,8,6,4,7,7,3,4,2,8,3,7,4,8,6,4,7,7,
4,5,3,7,4,6,5,7,5,3,6,6,2,3,1,7,2,6,3,7,7,5,6,6, 4,5,3,7,4,6,5,7,5,3,6,6,
3,6,4,8,5,7,4,8,7,7,7,7,3,4,2,8,3,7,4,8,6,4,7,7, 4,5,3,7,4,6,5,7,5,3,6,6,
6,7,5,7,6,6,7,7,5,3,6,6,6,7,5,7,6,6,7,7,5,3,6,6, 5,6,6,6,5,7,8,6,4,2,5,5,
6,7,5,7,6,6,7,7,5,3,6,6, 6,7,5,7,6,6,7,7,5,3,6,6, 3,4,4,8,3,9,6,8,6,4,7,7,
3,4,2,6,3,5,4,6,6,4,5,5, 4,5,3,5,4,4,5,5,5,3,4,4, 5,6,4,8,5,7,6,8,6,4,7,7,
5,6,4,8,5,7,6,8,6,4,7,7, 4,5,5,7,4,8,7,7,5,3,6,6, 5,6,4,8,5,7,6,8,6,4,7,7,
3,4,2,6,3,5,4,6,6,4,5,5, 5,6,4,8,5,7,6,8,6,4,7,7, 5,6,4,8,5,7,6,8,6,4,7,7,
5,4,4,6,5,5,6,6,6,6,7,5, 6,5,5,5,6,4,5,5,5,7,6,4, 4,5,3,7,4,6,5,7,7,5,6,6,
6,5,5,5,6,4,5,5,5,7,6,4, 3,4,4,6,3,7,6,6,6,4,5,5), nrow=12))
rownames(obj.vals) <- c("Enclosed Offices", "Open Office", "Reception / Greeter", "Studio / Classroom",
                        "Conference / Meeting Room", "Gallery", "Public / Lobby / Waiting", 
                        "Collaborative Space", "Mechanical / Support", "Storage / Archives", 
                        "Fabrication", "Performance")

For each of the 12 programs, let's set a maximum number of repetitions that can be assigned to all rooms combined. Note this is something that I added, since this data was not provided by the OP. (Restrict them from getting assigned to too many rooms.)

max_programs <- c(1,2,3,1,5,2,3,4,1,3,1,2)
    
library(lpSolveAPI)

nrooms <- 28
nprgs <- 12
ncol = nrooms*nprgs

lp_matching <- make.lp(ncol=ncol)
#we want integer assignments
set.type(lp_matching, columns=1:ncol, type = c("integer"))

# sum r,p Crp * Xrp
set.objfn(lp_matching, obj.vals) #28 rooms * 12 programs
lp.control(lp_matching,sense='max')

#' Set Max Programs constraints
#' No more than max number of programs over all the rooms
#' X1p + x2p + x3p ... + x28p <= max(p) for each p
Add_Max_program_constraint <- function (prog_index) {
  prog_cols <- (0:(nrooms-1))*nprgs + prog_index
  add.constraint(lp_matching, rep(1,nrooms), indices=prog_cols, rhs=max_programs[prog_index])
}
#Add a max_number constraint for each program
lapply(1:nprgs, Add_Max_program_constraint)

#' Sum of all the programs assigned to each room, over all programs 
#' area_1 * Xr1+ area 2* Xr2+ ... + area12* Xr12 <= room.size[r] for each room
Add_room_size_constraint <- function (room_index) {
  room_cols <- (room_index-1)*nprgs + (1:nprgs) #relevant columns for a given room
  add.constraint(lp_matching, xt=program.size, indices=room_cols, rhs=room.size[room_index])
}
#Add a max_number constraint for each program
lapply(1:nrooms, Add_room_size_constraint)

To solve this:

> solve(lp_matching)
> get.objective(lp_matching)
 [1] 195
get.variables(lp_matching) # to see which programs went to which rooms

> print(lp_matching)
Model name: 
  a linear program with 336 decision variables and 40 constraints

You can also write the IP model to a file to examine it:

#Give identifiable column and row names
rp<- t(outer(1:nrooms, 1:nprgs, paste, sep="_"))
rp_vec <- paste(abc, sep="")
colnames<- paste("x_",rp_vec, sep="")
# RowNames
rownames1 <- paste("MaxProg", 1:nprgs, sep="_")
rownames2 <- paste("Room", 1:nrooms, "AreaLimit", sep="_")
dimnames(lp_matching) <- list(c(rownames1, rownames2), colnames)

write.lp(lp_matching,filename="room_matching.lp")

Hope that helps.

Update based on follow-up questions

Follow up Question 1: Alter the code to ensure that every room has at least one program.

Add the following constraint set:

X_r_p >= 1 for all r

Note: Since this is a maximization problem, the optimal solution should honor this constraint by default. In other words, it will always assign a program to any room if it could, assuming positive scores for assigning.

Follow up question 2: Another question is if I can ask it to have more than 28 programs in total? For instance, if I want 28 Enclosed Offices they almost all could fit in the one room with 2938 square feet. How then can I ask R to still find other programs if the max is set for 28?

In order to achieve this goal, you can do it a bit differently. Do not have the sum of all programs <= 28 constraint at all. (If you note in the solution above, my constraints are slightly different.)

The constraint:

 Sum(r) X_r_p <= Max_allowable(p) for each program p

only limits the max per type of program. There is no limit on the total. Also, you don't have to write one such constraint for each type of program. Write this constraints only if you want to limit its occurrence.

To generalize this, you can set the lower and upper bounds for the total of each type of programs. This will give you very fine control over the assignments.

min_allowable(p) <= sum(over r) X_r_p <= max_allowable(p) for any program type p