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 usinggurobipy
, 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 orgurobipy
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%
gurobipy
takes 174.73 seconds andcvxpy
takes 54.93 seconds when I solved the ILP for 300 times (see the Question section). $\endgroup$