Skip to main content

Quadratic Objective Function with Linear Constraints Case

Answered
Muhamad Fikri
  • Gurobi-versary
  • Curious
  • Conversationalist

Hi,

I am new to this optimization topic and currently learning about optimization with quadratic objective functions.

In this case, I want to learn how to assign a limited number of people (with different sizes) to a room (limited size). I have an objective function

\[\min \sum_{j = 1}^m ( room\_size[j] - (\sum_{i=1}^n people\_size[i] x_{ij}))^2\]

where \(x\) is the decision variable. Given the constraints subject to

\(\sum_{i=1}^n x_{ij} >= 1\) for all \(j = 1, ..., m\) indicate I want to assign a room that can be filled by more than 1 person

\(\sum_{j=1}^m x_{ij} >= 1\) for all \(i = 1, ..., n\) indicate a person only enter 1 room.

with \(x_{ij}\) is Binary.

Then, I put the formulation into a code (below)

import numpy as np
import gurobipy as gp
from gurobipy import GRB

# Initial condition
people_size = np.array([18., 15., 11., 9., 14., 13.])
room_size = [60., 64.50, 54.50]

N = len(people_size) # number of items
M = len(room_size) # number of available rolls

assignment_matrix = np.zeros((N, M)) # put the assignment matrix
model = gp.Model("QP_General_Assignment_Problem")

# Create the decision variable the x_ij
x_ij = model.addVars(N, M, vtype = GRB.BINARY, name = 'x_ij')

# Create constraints
model.addConstrs((gp.quicksum(x_ij[i,j] for i in range(N)) >= 1) for j in range(M))
model.addConstrs((gp.quicksum(x_ij[i,j] for j in range(M)) == 1) for i in range(N))

# Create the objective function sum(room_size [j] - sum(people_size[i] x_ij) ** 2
model.setObjective(gp.quicksum(room_size [j] - gp.quicksum(people_size[i] * x_ij[i, j]  for i in range(N) ) for j in range(M))**2, sense = GRB.MINIMIZE)

# Find the optimal solution
model.optimize()
model.write("people_assignment.lp")

# Create the assignment matrix
assignment_matrix = np.zeros((N, M))
for i in range(N):
  for j in range(M):
    if x_ij[i, j].x == 1:
       assignment_matrix[i, j] = x_ij[i, j].x 
print(assignment_matrix)

 

From the code above, the result will be a matrix with 1 as an indication of people_i assigned to room_j. The code works, but I am not sure whether the quadratic in my objective function is correct or not. Somehow it is too good to be true. Therefore, could you suggest to me if the code I wrote above is a proper one to represent the quadratic objective function? Also, I believe the formula that I have is Quadratic Programming with Binary decision variables, but is the approach I use in the code correct?

Thank you for your valuable time in giving me insight

0

Comments

10 comments

Date Votes
  • Riley Clement
    • Gurobi Staff
    Edited

    Hi Muhamad,

    Well done on what you have put together.  The objective function is almost correct, you just need to bring the square inside the sum, so that you have a sum of squares, not a square of sums.  Auxiliary variables will help keep Gurobi happy, i.e.

    aux = model.addVars(M)
    for j in range(M):
        model.addConstr(aux[j] == room_size[j] - gp.quicksum(people_size[i] * x_ij[i, j]  for i in range(N) ))
    model.setObjective(gp.quicksum(aux[j]**2 for j in range(M)), sense = GRB.MINIMIZE)

    - Riley

    1
  • Muhamad Fikri
    • Gurobi-versary
    • Curious
    • Conversationalist

    Hi Riley

    Thank you, for your help in making my code correctly working.
    Previously I was confused when I tried to modify the room size number, but the assignment matrix was not changed at all.

     

    I do have another question. What is the meaning of this?

    Cutting planes:
      Gomory: 3
      MIR: 2
      Zero half: 9
      RLT: 5
      BQP: 10

     

    Thank you

     

    0
  • Riley Clement
    • Gurobi Staff

    Hi Muhamad,

    This section of the log summarizes which cutting planes were applied and how many.

    Please see the Cutting Planes section of MIP - A Primer on the Basics  and Cutting Plane Methods for Solving MIP Problems for more information.

    Some of our parameters, such as MIRCuts, give you some control over these cuts individually (as opposed to the broad parameter Cuts)

    - Riley

    0
  • Muhamad Fikri
    • Gurobi-versary
    • Curious
    • Conversationalist

    Hi Riley,

    The reference is good, I will read them all.

    Actually, I have a question again. Since I am already done with the Quadratic Objective function using normal equation. Then, I tried to extend the use of Gurobi using matrix form but I got a problem as it did not produce the same result as the previous one. Can I write the question here or should I create a new section?

     

    0
  • Riley Clement
    • Gurobi Staff

    Hi Muhamad,

    Seems like it is related to the original question, so feel free to write here and I'll try to get to it tomorrow.

    - Riley

    0
  • Muhamad Fikri
    Edited
    • Gurobi-versary
    • Curious
    • Conversationalist

    Hi, Riley

    Thank you for your help before.

    The idea I use a matrix is to have a better understanding of the equation form of my objective function. If my understanding is correct, the matrix form will be (let's assume room_size = r, and people_size =p) in vector and x_ij is also a vector.

    \[q = (r - (I_m \otimes p^T) x)^2\]

    \[q = (r - (I_m \otimes p^T) x)^T (r - (I_m \otimes p^T) x)\]

    \[q = r^T r - 2 r^T \tilde{p} x + x^T \tilde{p}^T \tilde{p} x\]

    Here \Tilde{p} = (I_m \otimes p^T), and r^Tr is constant we can ignore it. Please correct me if my calculation is incorrect. Then, the q will be my objective function. Since in QP matrix form, it has a form as

    \[1/2 x^TQx + c^T x\]

    Should I change the q objective function in my code? My current perspective is I need to split it and decide the value of Q and c^T, where Q = 2* x^T \Tilde{p}^T \Tilde{p} x, and c^T = 2 r^T \Tilde{p} x

    The code form is

    Q = 2 * np.kron(np.eye(N), people_size.T).T @ np.kron(np.eye(N), people_size.T)
    c_trans = -2 * people_size.T @ np.kron(np.eye(N), people_size.T)

    When I combined Q and c_trans as one variable let's say obj_func and input it to

    model.setObjective(obj_fun, GRB.MINIMIZE)

    I ended up with the error

       model.setObjective(obj_fun, GRB.MINIMIZE)
      File "src\gurobipy\model.pxi", line 1448, in gurobipy.Model.setObjective
      File "src\gurobipy\util.pxi", line 44, in gurobipy.__simpleexpr
    gurobipy.GurobiError: Unable to convert argument to an expression

    However, when I put them into the code as below, it worked but the result was not the same as the equation version result before.

    import numpy as np
    import gurobipy as gp
    from gurobipy import GRB

    # Define the people size and room size
    people_size = np.array([18., 15., 11., 9., 14., 13.])
    room_size = np.array([60., 64.50, 54.50])

    N = len(people_size) 
    M = len(room_size)

    res = 2 * np.kron(np.eye(M), people_size.T).T @ np.kron(np.eye(M), people_size.T)
    c_res = -2 * room_size.T @ np.kron(np.eye(M), people_size.T)

    # This part to handle Q for being not PSD or PD
    # Q' = Q - \lambda_min I

    rows, cols = res.shape
    eig = (np.linalg.eigvalsh(res))
    cvx_res = res - np.min(eig) * np.eye(rows)

    # The obj_func just for trial
    obj_fun = res + c_res

    # Create a Gurobi model
    model = gp.Model("QuadraticObjective")

    # Define decision variables
    x = {}
    for i in range(N):
        for j in range(M):
            x[i, j] = model.addVar(vtype=GRB.INTEGER, name=f"x_{i}_{j}")

    # Set objective function 0.5 x^T Q x
    obj = gp.QuadExpr()
    for i in range(N):
        for j in range(N):
            for k in range(M):
                for l in range(M):
                    obj += 0.5 * x[i, k] * cvx_res[i, j] * x[j, l]

    # Add c^T x to the objective function
    for i in range(N):
        for j in range(M):
            obj += 0.5 * x[i, j] * c_res[i]

    model.setObjective(obj_fun, GRB.MINIMIZE)


    # Add constraints here if needed
    model.addConstrs((gp.quicksum(x[i,j] for i in range(N)) >= 1) for j in range(M))
    model.addConstrs((gp.quicksum(x[i,j] for j in range(M)) == 1) for i in range(N))

    # Optimize the model
    model.optimize()

    assignment_matrix = np.zeros((N, M))
    for i in range(N):
      for j in range(M):
        if x[i, j].x == 1:
          assignment_matrix[i, j] = x[i, j].x * people_size[i]
    print(assignment_matrix)


    Sorry, I put a new Q to avoid the Q from not being PSD or PD. Do you have any suggestions on this matrix form? Although I believe you will have many insights I am not so sure whether I should put the constraints in matrix form or not.

     

    Thank you

    0
  • Riley Clement
    • Gurobi Staff
    Edited

    Hi Muhamad,

    A formulation with the matrix API will look like this:

    model = gp.Model("QuadraticObjective")
    x = model.addMVar((N, M), vtype="B")
    model.addConstr(x.sum(axis=0) >= 1)
    model.addConstr(x.sum(axis=1) == 1)
    model.setObjective((r - p@x)@(r - p@x))
    model.optimize()

    #show solution
    print(x.X)

    where p = people_size, and r = room_size

    Please let me know if you have any questions.  Note that in putting this together I realized that the objective function I proposed earlier was not correct, and the wrong part of the expression was being squared.  I have edited it to fix.

    - Riley

    0
  • Muhamad Fikri
    • Gurobi-versary
    • Curious
    • Conversationalist

    Hi Riley

    Thank you so much for your prompt reply.

    The code from you was really. I did not expect that I just needed to put it that way.

    My assumption from the matrix that I have is that it will have a nonconvex situation if I am not wrong. How does Gurobi tell us that our model is convex or nonconvex?

    0
  • Riley Clement
    • Gurobi Staff

    Hi Muhamad,

    If there are nonconvex constraints or objective then Gurobi will not let the optimization proceed unless the NonConvex parameter is set to 2 (and will issue an error message to this effect).

    You can find a bit more information in the Quadratic Constraints section in this Constraints section of our Reference Manual.

    - Riley

    0
  • Muhamad Fikri
    • Gurobi-versary
    • Curious
    • Conversationalist

    Hi Riley

    Thank you very much once again for your help.

    It clears now the problem that I have.

     

    0

Please sign in to leave a comment.

Didn't find what you were looking for?

New post