为什么 CVXOPT 会为这种非线性网络流优化给出等级错误?
Why does CVXOPT give a rank error for this nonlinear network flow optimisation?
我正在考虑使用 cvxopt 来解决一些非线性网络流量优化问题。为了了解基础知识,我尝试使用一个非常简单的测试网络,只有 4 个顶点和 5 个边。
我的网络看起来像 this。蓝色和红色节点分别是源和汇。
每条边的成本是:
alpha*x**2
其中 x 是包含每条边上的流的向量,alpha 是某个系数。那么我的优化问题是:
min sum(alpha*x**2)
subject to:
E*x = b
x >= 0
其中 E 是边弧入射矩阵,b 是包含源和汇的向量。因此,矩阵向量方程 Ex = b 应该只执行基尔霍夫定律(即每个节点处的流量守恒)。
我的 python 脚本是:
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']
我得到的错误是:
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
我不确定如何从这里开始。我假设这个优化问题可以用 cvxopt 解决,因为它很简单,可以手动找到最佳流程。如果有人能告诉我如何纠正这段代码,或者告诉我为什么这类问题不适合这个软件,我将不胜感激。
提前致谢。
仔细想了想,发现这个问题是因为cvxopt要求矩阵在等式约束中的秩不小于等式约束的个数。
在我的例子中,这意味着我的关联矩阵的等级必须等于网络中的节点数。然而,图论中的一个结果是任何具有 n 个节点的简单连通图都将具有秩为 n-1 的关联矩阵。这会产生排名错误。
我解决这个问题的方法是选择一个节点,并向其添加两条额外的边:一条从该节点开始但无处可去,另一条不知从何而来并终止于该节点。这实际上向矩阵添加了两列。然后我在矩阵中添加了一行,要求这两条新边上的流量之和为零。
通过这种方式,我增加了矩阵的秩,而没有添加任何额外的节点。原始网络上的流量不受添加这些边的影响,因为我要求新边上的流量保持为零。
这种方法有点笨拙,但似乎可以解决问题。
运行今天进入同样的问题。
另一种解决方法是删除冗余约束。
获取矩阵 E 的 SVD:
E = U S V'
S 的形状与 E 相同,最后一行全为零(因为您的矩阵是秩 3)。
让 y = V' x
并重新安排你的等式约束
E x = b
U S V' x = b
S V' x = U' b
最后一行U' b
必须为零,否则问题无法解决
设T
为S V'
的前3行,c
为U' b
的前3行。
那么就可以使用等式约束了
T x = c
或者,cvxpy
does a similar thing 使用 QR 分解
我正在考虑使用 cvxopt 来解决一些非线性网络流量优化问题。为了了解基础知识,我尝试使用一个非常简单的测试网络,只有 4 个顶点和 5 个边。
我的网络看起来像 this。蓝色和红色节点分别是源和汇。
每条边的成本是:
alpha*x**2
其中 x 是包含每条边上的流的向量,alpha 是某个系数。那么我的优化问题是:
min sum(alpha*x**2)
subject to:
E*x = b
x >= 0
其中 E 是边弧入射矩阵,b 是包含源和汇的向量。因此,矩阵向量方程 Ex = b 应该只执行基尔霍夫定律(即每个节点处的流量守恒)。
我的 python 脚本是:
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']
我得到的错误是:
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
我不确定如何从这里开始。我假设这个优化问题可以用 cvxopt 解决,因为它很简单,可以手动找到最佳流程。如果有人能告诉我如何纠正这段代码,或者告诉我为什么这类问题不适合这个软件,我将不胜感激。
提前致谢。
仔细想了想,发现这个问题是因为cvxopt要求矩阵在等式约束中的秩不小于等式约束的个数。
在我的例子中,这意味着我的关联矩阵的等级必须等于网络中的节点数。然而,图论中的一个结果是任何具有 n 个节点的简单连通图都将具有秩为 n-1 的关联矩阵。这会产生排名错误。
我解决这个问题的方法是选择一个节点,并向其添加两条额外的边:一条从该节点开始但无处可去,另一条不知从何而来并终止于该节点。这实际上向矩阵添加了两列。然后我在矩阵中添加了一行,要求这两条新边上的流量之和为零。
通过这种方式,我增加了矩阵的秩,而没有添加任何额外的节点。原始网络上的流量不受添加这些边的影响,因为我要求新边上的流量保持为零。
这种方法有点笨拙,但似乎可以解决问题。
运行今天进入同样的问题。
另一种解决方法是删除冗余约束。
获取矩阵 E 的 SVD:
E = U S V'
S 的形状与 E 相同,最后一行全为零(因为您的矩阵是秩 3)。
让 y = V' x
并重新安排你的等式约束
E x = b
U S V' x = b
S V' x = U' b
最后一行U' b
必须为零,否则问题无法解决
设T
为S V'
的前3行,c
为U' b
的前3行。
那么就可以使用等式约束了
T x = c
或者,cvxpy
does a similar thing 使用 QR 分解