Efficient way of adding large number of constraints
I am experiencing slowness using for loop within addConstrs method. The python code is as followed:
from gurobipy import Model, quicksum, GRB
import numpy as np
m = Model('sparse_decision_variables')
num_dim1 = 100
num_dim2 = 200
num_constraints = 1000
feasible_arcs = [(i, j) for i in range(num_dim1) for j in range(num_dim2) if i + j % 5 != 3]
A = np.random.rand(num_dim2, num_constraints)
rhs = np.random.rand(num_constraints, 1)
# decision variable x is a sparse 10 * 10 matrix (addMVar doesnt support adding variables from list)
x = m.addVars(feasible_arcs, vtype=GRB.BINARY, name="x")
m.addConstrs((((quicksum(x[i, j] * A[j, k] for i, j in feasible_arcs)) <= rhs[k]) for k in range(num_constraints)), name='slow_constraint' )
m.update()
m.optimize()
Suppose a and b are the number of rows and columns of decision variable x, c is the size of rhs (number of constraints). For clearity, we use @ to represent matrix multiplication and * to represent element wise multiplication.
In this example, a = 100, b=200, c=1000
In matrix form, x @ A has dimension a*c (x is a a * b matrix and A is a b * c matrix). The constraint basically wants to express row sum of x@ A <= rhs. MATLAB allows to add constraints like one vector containing decision variables is less than another vector. But in python gurobi API, all I know is to loop over each element in the vector (for k in range(num_constraints)).
We have already tried:
m.addConstrs((((quicksum(x[i, j] * A[j, :] for i, j in feasible_arcs)) >= rhs)), name='slow_constraint' )
m.addConstrs((((np.dot(x[i, j] , A[j, :] for i, j in feasible_arcs)) >= rhs)), name='slow_constraint' )
But neither worked.
Could anyone propose an idea to efficiently add this type of constraints? Thanks

We have already tried:
m.addConstrs((((quicksum(x[i, j] * A[j, :] for i, j in feasible_arcs)) >= rhs)), name='slow_constraint' )
m.addConstrs(((((np.dot(x[i, j] , A[j, :]) for i, j in feasible_arcs)) >= rhs)), name='slow_constraint' )
I should have added a parenthesis closing np.dot in the second constraint but it does not generate a valid constraint as we expected

Hi Shiyao,
If I understand correctly, you want to add the constraint \( XAe \leq d \), where \( X \in \mathbb{R}^{m \times n} \) is a matrix of variables, \( A \in \mathbb{R}^{n \times p} \) is a coefficient matrix, \( e \in \mathbb{R}^p \) is the vector of ones in dimension \( p \), and \( d \in \mathbb{R}^p \) is your RHS vector. If so, the rowwise sum of \( XA \) (i.e., \( XAe \)) is an \( m \)dimensional vector, which doesn't match the RHS dimension of \( p \). Could you clarify which constraints you're trying to add? Maybe this should be the columnwise sum of \( X A \)?
Also, is the following condition correct when creating the feasible_arcs list?
i + j % 5 != 3
I wonder if there should be parentheses around i + j. As written, the feasible_arcs list is not very sparse, containing 19840 out of a possible 20000 arcs.
Thanks,
Eli

Hi Eli,
Thank you so much for you reply!
I do appologize for typos in code and type settings. We are trying to find out how to type latex here.
For your first question regarding to the dimension mismatch, you are absolutely correct, I mean column sum. Following your notation, I tried to add a constraint e^T XA <= d, where e is a n * 1 column vector of ones.The dimension of X, A and d are exactly as you stated above.
For the feasible_arcs list, you are right again. It should be parenthesised like
(i + j) % 5 != 3
For decision variables, I could have used addMVar method to creat a two dimensional decision variable. But in our actual problem, len(feasible_arcs) is less than 1% of m*n. So we want to create decision variables from list to reduce problem size.
The reason why I did not create very sparse decision variables is that I want to show that adding constraint takes a significant amount of time. Alternatively, we could increase the size of X and A and in the list,
(i + j) % 100 == 1
Shiyao

Hi Shiyao,
Thanks for the clarification. I think you could add these constraints more efficiently if you had a mapping that defines the inbound arcs to a given node. E.g.:
inbound = {}
for i,j in feasible_arcs:
if j in inbound:
inbound[j].append(i)
else:
inbound[j] = [i]Then, you can add the constraints as follows:
m.addConstrs((quicksum(A[j,k] * quicksum(x[i,j] for i in inbound[j]) for j in inbound) <= rhs[k]) for k in range(num_constraints))
Note that because the matrix \( A \) is completely dense in your code, the number of nonzeros in the coefficient matrix is equal to the number feasible arcs times the number of constraints. Thus, with the full 20000 arcs, the coefficient matrix contains 20 million nonzeros. Hopefully the model building is faster for more sparse arc sets.
Eli
Please sign in to leave a comment.
Comments
4 comments