I'm trying to find a library for convex optimization in python which can handle a large amount of variables and constraints.
My problem is in the form
$$ \min p^T \cdot q + |c| \\ \text{lower bound} \leq q \leq \text{upper bound} \\ \text{lower bound} \leq c \leq \text{upper bound} \\ Aq + Bc = d \\ \text{lower bound} \leq Rq \leq \text{upper bound} $$
Where $p, d$ are vectors and $A, B$ are sparse matrices and $q, c$ are variables. $R$ is a special matrix which might be the cause of my issue.
$$ Rq = \\ \begin{pmatrix} q_1 \\ q_2 \\ q_3 \\ q_4 \\ q_1 + q_5 \\ q_2 + q_6 \\ q_3 + q_7\\ q_4 + q_8 \\ q_1 + q_5 +q_9 \\ q_2 + q_6 + q_{10} \\ q_3 + q_7 + q_{11} \\ q_4 + q_8 + q_{12} \\ \vdots \\ \sum_{k=0}^n q_{4k+1} \\ \sum_{k=0}^n q_{4k+2} \\ \sum_{k=0}^n q_{4k+3} \\ \sum_{k=0}^n q_{4k+4} \end{pmatrix} \\ =\\ \begin{pmatrix} 1 & 0 & \dots & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\ 0 & 1 & 0 & \dots & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\ 0 & 0 & 1 & 0 & \dots & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\ 0 & 0 & 0 & 1 & 0 & \dots & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\ 1 & 0 & 0 & 0 & 1 & 0 & \dots & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\ 0 & 1 & 0 & 0 & 0 & 1 & 0 & \dots & 0 & 0 & 0 & 0 & 0 & 0 & 0\\ 0 & 0 & 1 & 0 & 0 & 0 & 1 & 0 & \dots & 0 & 0 & 0 & 0 & 0 & 0\\ 0 & 0 & 0 & 1 & 0 & 0 & 0 & 1 & 0 & \dots & 0 & 0 & 0 & 0 & 0\\ 1 & 0 & 0 & 0 & 1 & 0 & 0 & 0 & 1 & 0 & \dots & 0 & 0 & 0 & 0\\ 0 & 1 & 0 & 0 & 0 & 1 & 0 & 0 & 0 & 1 & 0 & \dots & 0 & 0 & 0\\ 0 & 0 & 1 & 0 & 0 & 0 & 1 & 0 & 0 & 0 & 1 & 0 & \dots & 0 & 0\\ 0 & 0 & 0 & 1 & 0 & 0 & 0 & 1 & 0 & 0 & 0 & 1 & 0 & \dots & 0\\ \vdots & \vdots & \vdots & \vdots & \vdots & \vdots & \vdots & \vdots & \vdots & \vdots & \vdots & \vdots & \vdots & \vdots & \vdots & \end{pmatrix}q $$
code
from cvxopt.modeling import variable, op, sum
q = variable(numb_prod, 'q')
c = variable(numb_capacity, 'c')
prod_price = price_purchase_power.trans() * q
problem = op(sum(prod_price) + 1e-6*sum(abs(c)), [constraints])
problem.solve(format='sparse', solver='glpk', options=
{'glpk'{'msg_lev':'GLP_MSG_OFF'}})
error
Traceback (most recent call last): File "C:\Users\OsRas\AppData\Local\Continuum\anaconda3\lib\site-packages\IPython\core\interactiveshell.py", line 3325, in run_code exec(code_obj, self.user_global_ns, self.user_ns) File "", line 489, in q, c, total_cost_x, water_values, list_numb_of_capacity_variables = solution(demand, price_purchase_power, supply) File "", line 391, in solution test = problem.solve(format='sparse', solver='glpk', options={'glpk':{'msg_lev':'GLP_MSG_OFF'}}) File "C:\Users\OsRas\AppData\Local\Continuum\anaconda3\lib\site-packages\cvxopt\modeling.py", line 2600, in solve t = self._inmatrixform(format) File "C:\Users\OsRas\AppData\Local\Continuum\anaconda3\lib\site-packages\cvxopt\modeling.py", line 2479, in _inmatrixform for v in variables: vmap[v] = x[vslc[v]] File "C:\Users\OsRas\AppData\Local\Continuum\anaconda3\lib\site-packages\cvxopt\modeling.py", line 240, in getitem return (+self).getitem(key) File "C:\Users\OsRas\AppData\Local\Continuum\anaconda3\lib\site-packages\cvxopt\modeling.py", line 931, in getitem f._linear = self._linear[l] File "C:\Users\OsRas\AppData\Local\Continuum\anaconda3\lib\site-packages\cvxopt\modeling.py", line 1385, in getitem for v,c in iter(self._coeff.items()): OverflowError: Python int too large to convert to C long
Currently I'm having 108575 variables which is too large for cvxopt. I was hopping running an optimization which can handle the amount of variables in the millions...
Is there any library that can handle that? Or does one need to look into another language?
I can sacrifice some accuracy and make this problem LP by removing $|c|$ if that opens up other possibilities.
I've asked this question on stackoverflow as well. But perhaps there is a better auidience here.
Edit As the problem is formulated now we have
$\dim(A) = 26058 \times 73831 $ with average three $1$:s per row, rest zero, upper triangle.
$\dim(B) =26058 \times 34744$ with average of three elements per row, $1$:s and $-1$:s, rest zero, elements around diagonal and some one the diagonal.
$\dim(R) = 17280 \times 73831$
On Stackoverflow you shared that $|c|$ means $||c||_1$, so the problem can be made linear. I can think of two solution methods.
The first one is to use an interior point algorithm. Due to the size of the problem, the location of the $1$s is important. If either $AA^T + BB^T$ or $$\begin{pmatrix}A^TA & O \\ O & B^TB\end{pmatrix}$$ has a sparse Cholesky factorization (for example, it can be row/column permuted to an arrow matrix), this approach can be memory efficient and fast. You could try CPLEX or Gurobi. They should tell you the number of nonzeros in the cholesky factor.
The second option is to use ADMM. You could solve the subproblems in cvxpy or use a specialized package. I have never had much luck with ADMM, but your mileage may vary.