# Background

I am trying to implement the reranking algorithm proposed in the paper (Eqn 3). The algorithm is cast into an integer linear program which is trying to find permutation matrix $$\mathbf{X} \in \mathbb{R}^{n\times n}$$, where $$\mathbf{X}_{ij} =1$$ means ranking item $$i$$ to position $$j$$.

# Question

Initially, I implemented this algorithm using cvxpy and setting the solver to be Gurobi. The runtime is not satisfying. So I decided to directly call Gurobi Python API (i.e. gurobipy) rather than using the modeling interface (i.e. cvxpy). I expected the runtime to be improved as I skip the abstraction from cvxpy. However, the runtime is surprisingly much worse (see below).

# gurobipy runtime in seconds
174.727278
# cvxpy runtime in seconds
54.925721


The only two explanations I could think are

• My implementation based on gurobipy is not efficient.

• cvxpy somehow optimizes the formulation by using some transformation. This transformation makes it easier for Gurobi to discover the structure in the problem. However, when directly code the formulation using gurobipy, this structure is not well exploited by Gurobi.

My questions are

• Does any of my explanations make sense? Why this odd behavior might happen.
• What could I do to optimize my cvxpy code or gurobipy code to improve the runtime?

# Code to reproduce my result

import numpy as np
import cvxpy as cp
import gurobipy as gp

from datetime import datetime
from gurobipy import GRB

def runGUROBIImpl(n, theta=0.5, nIter=10):
for i in range(nIter):
with gp.Env(empty=True) as env:
env.setParam("OutputFlag", 0)
env.start()
with gp.Model(env=env) as model:
C = np.random.rand(n, n)
DCG = np.random.rand(n, n)
IDCG = np.random.rand()

P = model.addVars(n, n, vtype=GRB.BINARY)
obj = sum([C[i, j] * P[i, j] for i in range(n) for j in range(n)])
constrLHS = sum([DCG[i, j] * P[i, j] for i in range(n) for j in range(n)])

model.setObjective(obj, GRB.MINIMIZE)
model.addConstr(constrLHS >= theta * IDCG)
model.addConstrs((P.sum("*", i) == 1 for i in range(n)))
model.addConstrs((P.sum(i, "*") == 1 for i in range(n)))

model.optimize()

if not (model.status == GRB.OPTIMAL): print("unsuccessful...")

def runCVXPYImpl(n, theta=0.5, nIter=10):
C, DCG, IDCG = cp.Parameter((n, n)), cp.Parameter((n, n)), cp.Parameter()
P = cp.Variable((n, n), boolean=True)

objective = cp.Minimize(cp.sum(cp.multiply(C, P)))
constraints = [cp.sum(cp.multiply(DCG, P)) >= theta * IDCG,
cp.matmul(np.ones((1, n)), P) == np.ones((1, n)),
cp.matmul(P, np.ones((n, 1))) == np.ones((n, 1))]
problem = cp.Problem(objective, constraints)

for i in range(nIter):
C.value = np.random.rand(n, n)
DCG.value = np.random.rand(n, n)
IDCG.value = np.random.rand()
# solve problem
problem.solve(solver=cp.GUROBI, verbose=False)

if not (problem.status == cp.OPTIMAL): print("unsuccessful...")

np.random.seed(0)
start = datetime.now()
runGUROBIImpl(n=100, theta=0.5, nIter=300)
end = datetime.now()
print((end - start).total_seconds())

np.random.seed(0)
start = datetime.now()
runCVXPYImpl(n=100, theta=0.5, nIter=300)
end = datetime.now()
print((end - start).total_seconds())


# Logs

The algorithm proposed in the paper requires the linear program to be solved several hundreds times with updated parameter matrix each time to get the final result $$\mathbf{X}$$.

As it would be much too verbose to provide logs for several hundreds solving attempts, the following are the logs when I set nIter=1 for both runGUROBIImpl() and runCVXPYImpl().

Note when nIter=1 the situation for runtime is different from when nIter=300 as the cvxpy requires a lot of time to construct problem.

# gurobipy runtime in seconds
0.61062
# cvxpy runtime in seconds
15.098772


## runGUROBIImpl() Logs

Gurobi Optimizer version 9.1.1 build v9.1.1rc0 (linux64)
Thread count: 2 physical cores, 4 logical processors, using up to 4 threads
Optimize a model with 201 rows, 10000 columns and 30000 nonzeros
Model fingerprint: 0x4962970c
Variable types: 0 continuous, 10000 integer (10000 binary)
Coefficient statistics:
Matrix range     [2e-04, 1e+00]
Objective range  [7e-05, 1e+00]
Bounds range     [1e+00, 1e+00]
RHS range        [2e-01, 1e+00]
Found heuristic solution: objective 48.7825807
Presolve removed 1 rows and 0 columns
Presolve time: 0.03s
Presolved: 200 rows, 10000 columns, 20000 nonzeros
Variable types: 0 continuous, 10000 integer (10000 binary)

Root relaxation: objective 1.415358e+00, 210 iterations, 0.00 seconds

Nodes    |    Current Node    |     Objective Bounds      |     Work
Expl Unexpl |  Obj  Depth IntInf | Incumbent    BestBd   Gap | It/Node Time

*    0     0               0       1.4153579    1.41536  0.00%     -    0s

Explored 0 nodes (210 simplex iterations) in 0.03 seconds
Thread count was 4 (of 4 available processors)

Solution count 2: 1.41536 48.7826

Optimal solution found (tolerance 1.00e-04)
Best objective 1.415357907906e+00, best bound 1.415357907906e+00, gap 0.0000%


## runCVXPYImpl() Logs

Parameter OutputFlag unchanged
Value: 1  Min: 0  Max: 1  Default: 1
Changed value of parameter QCPDual to 1
Prev: 0  Min: 0  Max: 1  Default: 0
Gurobi Optimizer version 9.1.1 build v9.1.1rc0 (linux64)
Thread count: 2 physical cores, 4 logical processors, using up to 4 threads
Optimize a model with 201 rows, 10000 columns and 30000 nonzeros
Model fingerprint: 0x8b5dd50a
Variable types: 0 continuous, 10000 integer (10000 binary)
Coefficient statistics:
Matrix range     [2e-04, 1e+00]
Objective range  [7e-05, 1e+00]
Bounds range     [1e+00, 1e+00]
RHS range        [2e-01, 1e+00]
Found heuristic solution: objective 47.5478288
Presolve removed 1 rows and 0 columns
Presolve time: 0.03s
Presolved: 200 rows, 10000 columns, 20000 nonzeros
Variable types: 0 continuous, 10000 integer (10000 binary)

Root relaxation: objective 1.415358e+00, 210 iterations, 0.00 seconds

Nodes    |    Current Node    |     Objective Bounds      |     Work
Expl Unexpl |  Obj  Depth IntInf | Incumbent    BestBd   Gap | It/Node Time

*    0     0               0       1.4153579    1.41536  0.00%     -    0s

Explored 0 nodes (210 simplex iterations) in 0.03 seconds
Thread count was 4 (of 4 available processors)

Solution count 2: 1.41536 47.5478

Optimal solution found (tolerance 1.00e-04)
Best objective 1.415357907906e+00, best bound 1.415357907906e+00, gap 0.0000%

• Can you please show the logs? – RobPratt Jan 1 at 23:11
• @RobPratt I have provided the logs in the question. Please let me know if you need any additional info. – Mr.Robot Jan 1 at 23:47
• Both logs show the problem being solved almost immediately (0.03 seconds). Can you determine whether there is any significant difference between the cvxpy and gurobipy runs in the amount of time that they spend in the Gurobi solver engine? If not, then the possibility of an inefficient gurobipy implementation should be considered. Like cvxpy, gurobipy is a python-based modeling language, which can encounter inefficiencies if not written carefully. Concerning efficiency of various ways of building linear expressions in gurobipy, see gurobi.com/documentation/9.1/refman/py_lex.html. – 4er Jan 4 at 16:59
• @4er Thanks for pointing out the 0.03 second solving time. I did not notice this. And yes, the difference of runtime of two versions is substantial, where gurobipy takes 174.73 seconds and cvxpy takes 54.93 seconds when I solved the ILP for 300 times (see the Question section). – Mr.Robot Jan 5 at 5:04

## 2 Answers

As others have pointed out already, you are not solving the same instances. When writing out the MPS files using a gurobi.env file containing GURO_PAR_DUMP=1, we can see that the instances differ (here n=2, nIters=1; left is cvxpy, right is Gurobi):

To get the signs in order, you could change this line in the Gurobi method

model.addConstr(constrLHS >= theta * IDCG)


to

model.addConstr(-constrLHS <= -theta * IDCG)


Apparently, cvxpy changes this internally. This still leaves the variable order different in both models. I didn't manage to get this identical for both codes. I did verify that the only difference is in the ordering - the objective values coincide.

That Gurobi takes more time in some configurations (especially when the instance size is larger) is probably due to the fact that you are constructing a new model every time (you should reuse the environment, by the way). In cvxpy on the other hand, the model structure is only constructed once and then reused for all instances. This also doesn't require creating a new set of variables every time. I am not sure how this is handled internally, though, since at some point the model data needs to be handed over to Gurobi.

Note that the model fingerprint differs. I suspect that the variable or constraint orders are different.

• Thank you! I did not know the usage of model fingerprint. So in my case, the runGUROBIImpl() and runCVXPYImpl() are actually solving different problems and the latter is easier for Gurobi to solve? – Mr.Robot Jan 2 at 0:57
• Probably the same problem but with the variables or constraints in different orders. See also or.stackexchange.com/questions/4042/… – RobPratt Jan 2 at 1:09