An answer to this question on the Scientific Computing Stack Exchange.

Question

I'm trying to optimize a binary portfolio vector to be greater than a benchmark using CVXPY.

import cvxpy as cp
import numpy as np
# Generate a random non-trivial quadratic program.
n = 10 # number of options
np.random.seed(1)
mu = np.random.randn(n) # expected means
var_covar = np.random.randn(n,n) # variance-covariance matrix
var_covar = var_covar.T.dot(var_covar) # cont'd
bench_cov = np.random.randn(n) # n-length vector of cov(benchmark, returns)
lamd = 0.01 # risk tolerance
# Define and solve the CVXPY problem.
x = cp.Variable(n, boolean=True)
prob = cp.Problem(cp.Maximize(mu.T@x + lamd * (cp.quad_form(x, var_covar) - (2 * bench_cov.T@x))), [cp.sum(x) == 4])
prob.solve()

I get this error using CVXPY version 1.1.0a0 (downloaded directly from github):

DCPError: Problem does not follow DCP rules. Specifically:

The objective is not DCP, even though each sub-expression is.

You are trying to maximize a function that is convex.

From what I've read maximizing a convex function is very difficult, but I got this equation from a paper that solves the same equation using Gurobi's BQP solver. I figure I must be doing something wrong as I'm new to quadratic programming and CVXPY.

Thank you!

Answer

cvxpy's rules for disciplined convex programming are listed here.

Notably, it states that:

The DCP rules require that the problem objective have one of two forms:

Minimize(convex)

Maximize(concave)

Indeed, your program runs if we remove the quad_form and doesn't run if we leave only the quad_form:

prob = cp.Problem(cp.Maximize(cp.quad_form(x, var_covar)), [cp.sum(x) == 4])

So it doesn't work because you are not obeying the rules of disciplined convex programming.

How can you solve this problem?

1. Fix the math. Maybe you made a mistake?
2. cvxpy's implementation of DCP is built from atoms, so any problem you want to solve must be expressible in this atoms. If you believe your program should be solvable by the techniques used by cvxpy, then perhaps you can rejigger your math to use a different set of atoms to express the same problem. (I think this is unlikely given the simplicity of the problem.) This might also be why Gurobi can solve the problem: it may have access to a more expansive set of inter-atom manipulations.
3. Don't use disciplined convex programming. That boolean=True means that you're already living in integer programming land. To my knowledge, DCP can't help you here.

Edit:

Perhaps the matrix in the paper you are reading is negative definite. In that case, you'd be maximizing a concave function, which is allowed. In this case, you've simply built your example wrong. Try this instead:

import cvxpy as cp
import numpy as np
# Generate a random non-trivial quadratic program.
n = 10 # number of options
np.random.seed(1)
mu = np.random.randn(n) # expected means
# var_covar = np.random.randn(n,n) # variance-covariance matrix
# var_covar = var_covar.T.dot(var_covar) # cont'd
bench_cov = np.random.randn(n) # n-length vector of cov(benchmark, returns)
#Make a negative definite matrix
var_covar = -np.array(range(n))
var_covar = np.diag(var_covar)
lamd = 0.01 # risk tolerance
# Define and solve the CVXPY problem.
x = cp.Variable(n)
#Ensure matrix is negative definite so we are maximizing a concave function
assert np.all(np.linalg.eigvals(var_covar) <= 0)
prob = cp.Problem(cp.Maximize(cp.quad_form(x, var_covar)), [cp.sum(x) == 4])
prob.solve()