I'm considering using cvxopt to solve some nonlinear network flow optimisation problems. In order to understand the basics, I am trying it with a very simple test network with just 4 vertices and 5 edges.
My network looks like this. The blue and red nodes are a source and sink respectively.
The cost on each edge is:
alpha*x**2
where x is the vector containing the flow on each edge, and alpha is some coefficient. My optimisation problem is then:
min sum(alpha*x**2)
subject to:
E*x = b
x >= 0
where E is the edge-arc incidence matrix and b is the vector containing sources and sinks. The matrix-vector equation Ex = b should therefore just enforce Kirchoff's laws (i.e. conservation of flow at each node).
My python script to do this is:
from cvxopt import matrix, spdiag, solvers
#Four vertices, five edges
E = matrix(0.0, (4,5))
E[0,0] = 1.0
E[0,1] = 1.0
E[1,0] = -1.0
E[1,2] = 1.0
E[1,3] = 1.0
E[2,1] = -1.0
E[2,2] = -1.0
E[2,4] = 1.0
E[3,3] = -1.0
E[3,4] = -1.0
#Source-sink vector
b = matrix(0.0, (4,1))
b[0] = 0.5
b[3] = -0.5
#Coefficient in the quadtratic
alpha = 1.0
#Function to pass to cvxopt
def F(x=None, z=None):
if x is None:
return 0, matrix(0.05, (5,1))
if min(x) <= 0.0:
return None
f = sum(alpha*(x**2))
Df = (2.0*alpha*x).T
if z is None:
return f, Df
D2f = 2.0*alpha*matrix(1.0, (5,1))
H = spdiag(z[0]*D2f)
return f, Df, H
#Solve
x = solvers.cp(F, A=E, b=b)['x']
The error I get is:
pcost dcost gap pres dres
0: 0.0000e+00 1.2500e-02 1e+00 1e+00 2e-01
Traceback (most recent call last):
File "simple_network.py", line 45, in <module>
x = solvers.cp(F, A=E, b=b)['x']
File "/usr/local/lib/python2.7/dist-packages/cvxopt/cvxprog.py", line 1966, in cp
xdot_e, xaxpy_e, xscal_e, options = options)
File "/usr/local/lib/python2.7/dist-packages/cvxopt/cvxprog.py", line 782, in cpl
raise ValueError("Rank(A) < p or "\
ValueError: Rank(A) < p or Rank([H(x); A; Df(x); G]) < n
I am unsure how to proceed from here. I assumed this optimisation problem would be soluble with cvxopt, since it is simple enough to find the optimal flow by hand. I would appreciate it if somebody could tell me how to rectify this code, or else tell me why this type of problem is not suitable for this software.
Thanks in advance.
After thinking about it some more, I realise that this problem is caused because cvxopt requires that the rank of the matrix in the equality constraint be no less than the number of equality constrains.
In my case this means than the rank of my incidence matrix would have to be equal to the number of nodes in the network. However, a result in graph theory is that any simple, connected graph with n nodes will have an incidence matrix with rank n-1. This produces the rank error.
The way I got around this problem was to pick a node, and add two additional edges to it: one that starts at the node but leads nowhere, and another that comes from nowhere and terminates at the node. This in effect adds two columns to the matrix. I then added an additional row to the matrix to demand that the sum of the flow on these two new edges be zero.
In this way, I increase the rank of the matrix, without adding any additional nodes. The flow on the original network is not effected by the addition of these edges, since I demand that the flow on the new edges remains zero.
This is a bit of a hacky way to do it, but it seems to do the trick.