线性规划:我可以制定一个 objective 来同时最大化多个变量吗?

Linear Programming: can I formulate an objective to maximize multiple variables at once?

假设我有以下系统说明的一些变量和约束: 灰线可以拉伸和收缩由它们顶部的范围给定的量。蓝线只是端点,显示了灰线如何相互作用。

我的目标:我想使用线性规划来均匀地最大化灰色线条的大小,如图所示。你可以想象上面带有弹簧的灰色线条,它们都同样向外推。一个糟糕的解决方案是将所有蓝线尽可能地推到一边。请注意,此描述中有一点回旋余地,并且可能有多种解决方案 - 我所需要的只是让它们合理均匀,而不是让一个值最大化压扁其他所有值。

我试过的objective函数简单地最大化了线的总和:

maximize: (B - A) + (C - B) + (C - A) + (D - C) + (E - B) + (E - D) + (F - E) + (F - D) + (F - A) 

我很清楚这不是一个好的解决方案,因为项相互抵消并且一行中的增加只会使另一行中的它减少相同的量,所以 objective 永远不会加权在变量之间均匀分布最大化。

我还尝试最小化每条线与其中间可能范围的距离。对于行 B - A,其 (1,3) 范围内的中间值是 2。这是第一项的 objective:

minimize: |(B - A) - 2| + ...

为了实现绝对值,我将术语替换为 U 并添加了额外的约束:

minimize: U + ...
with: U <= (B - A - 2)
      U <= -(B - A - 2)

这与另一个问题相同objective:差异总是与另一行差异的变化成正比。我想如果我能对差进行平方就可以了,但我不能在线性求解器中输入它。

是否有一些 objective 函数可以实现我正在寻找的东西,或者线性求解器不是正确的工具?

我正在使用 Google OR-Tools,如果有帮助的话。

这里写出约束条件:

 1 <= B - A <= 3
 0 <= C - B <= 1
 1 <= C - A <= 4
 9 <= E - B <= 11
 7 <= D - C <= 11
 0 <= E - D <= 1
 3 <= F - E <= 5
 3 <= F - D <= 6
15 <= F - A < = 15

请记住,您最大的问题是您不知道自己到底想要什么。所以我不得不猜测。有时看到一些猜测可以帮助您完善您想要的内容,所以这对您来说并不算太糟糕,但它确实会使您的问题对于本网站的格式来说更加困难。

首先,我假设 spring 可以建模为有向无环图。也就是说,我可以将所有 spring 替换为指向右侧的箭头。永远不会有从右指向左的箭头(否则你的 spring 会弯曲成一个圆圈)。

完成后,您可以使用设置逻辑来确定最左边的蓝色条的标识。 (我假设只有一个 - 它留作练习以弄清楚如何概括。)然后您可以将此栏锚定在合适的位置。所有其他栏将相对于它定位。这个约束看起来像:

S[leftmost] = 0

现在,我们需要一些约束。

每条边i都有一个源点和一个终点(因为边是有向的)。调用源点位置S和终点位置E。此外,边缘具有最小长度l和最大长度L。由于我们固定了最左边的蓝色条的位置,因此连接到它的 spring 定义了它们的端点所在的间隔。这些端点是其他 springs, &c 的源点。因此,每条边对其端点的位置定义了两个约束。

S[i]+l[i] <= E[i]
E[i]      <= S+L[i]

顺便说一句,请注意我们现在可以制定一个简单的线性规划:

min 1
s.t.  S[leftmost]  = 0
      S[i]+l[i]   <= E[i]
      E[i]        <= S+L[i]

如果这个程序可以解决,那么你的问题就有了可行的解决方案。也就是说,条的长度不会对蓝条的位置产生相互不一致的描述。

现在,我们想要“均匀地最大化灰线的大小”,不管这是什么意思。

最小化与平均长度的偏差

这是一个想法。程序为每个小节选择的长度由 E[i]-S[i] 给出。让我们指定此长度应“接近”柱的平均长度 (L[i]+l[i])/2。因此,我们要为每个柱最小化的目标数量是:

(E[i]-S[i])-(L[i]+l[i])/2

有问题的是,这个值可以是正的也可以是负的,这取决于是否 (E[i]-S[i])>(L[i]+l[i])/2。这不好,因为我们想最小化与 (L[i]+l[i])/2 的偏差,该值应始终为正值。

为了解决这个问题,让我们对值进行平方,然后取平方根,得到:

sqrt(((E[i]-S[i])-(L[i]+l[i])/2)^2)

这似乎无法解决,但请耐心等待。

注意前面和取单元素向量的L2范数是一样的,所以我们可以改写为:

|(E[i]-S[i])-(L[i]+l[i])/2|_2

我们现在可以对每个柱的偏差求和:

|(E[0]-S[0])-(L[0]+l[0])/2|_2 + |(E[1]-S[1])-(L[1]+l[1])/2|_2 + ...

这给了我们以下优化问题:

min |(E[0]-S[0])-(L[0]+l[0])/2|_2 + |(E[1]-S[1])-(L[1]+l[1])/2|_2 + ...
s.t.  S[leftmost]  = 0
      S[i]+l[i]   <= E[i]
      E[i]        <= S+L[i]

这个问题用上面的形式不容易解决,但是我们可以通过引入一个变量来进行简单的操作t

min   t[0] + t[1] + ...
s.t.  S[leftmost]  = 0
      S[i]+l[i]   <= E[i]
      E[i]        <= S+L[i]
      |(E[i]-S[i])-(L[i]+l[i])/2|_2<=t[i]

这个问题与上一个问题完全一样。那么我们收获了什么?

优化是一种将问题转化为标准形式的游戏。一旦我们遇到标准形式的问题,我们就可以站在巨人的肩膀上,使用强大的工具来解决我们的问题。

上述操作已将问题变成 second-order cone problem (SOCP)。一旦变成这种形式,几乎可以直接解决。

这样做的代码如下所示:

#!/usr/bin/env python3

import cvxpy as cp
import networkx as nx
import matplotlib.pyplot as plt

def FindTerminalPoints(springs):
  starts = set([x[0] for x in springs.edges()])
  ends   = set([x[1] for x in springs.edges()])
  return list(starts-ends), list(ends-starts)

springs = nx.DiGraph()
springs.add_edge('a', 'b', minlen= 1, maxlen= 3)
springs.add_edge('a', 'c', minlen= 1, maxlen= 4)
springs.add_edge('a', 'f', minlen=15, maxlen=15)
springs.add_edge('b', 'c', minlen= 0, maxlen= 1)
springs.add_edge('b', 'e', minlen= 9, maxlen=11)
springs.add_edge('c', 'd', minlen= 7, maxlen=11)
springs.add_edge('d', 'e', minlen= 0, maxlen= 1)
springs.add_edge('d', 'f', minlen= 3, maxlen= 6)
springs.add_edge('e', 'f', minlen= 3, maxlen= 5)

if not nx.is_directed_acyclic_graph(springs):
  raise Exception("Springs must be a directed acyclic graph!")

starts, ends = FindTerminalPoints(springs)
if len(starts)!=1:
  raise Exception("One unique start is needed!")

if len(ends)!=1:
  raise Exception("One unique end is needed!")  

start = starts[0]
end   = ends[0]

#At this point we have what is essentially a directed acyclic graph beginning at
#`start` and ending at `end`

#Generate a variable for the position of each blue bar
bluevars = {n: cp.Variable(name=n) for n in springs.nodes()}
dvars    = {e: cp.Variable()       for e in springs.edges()}
#Anchor the leftmost blue bar to prevent pathological solutions
cons   = [bluevars[start]==0]
for s,e in springs.edges():
  print("Loading edge {0}-{1}".format(s,e))
  sv   = bluevars[s]
  ev   = bluevars[e]
  edge = springs[s][e]
  cons += [sv+edge['minlen']<=ev]
  cons += [ev<=sv+edge['maxlen']]
  cons += [cp.norm((ev-sv)-(edge['maxlen']-edge['minlen'])/2,2)<=dvars[(s,e)]]

obj  = cp.Minimize(cp.sum(list(dvars.values())))
prob = cp.Problem(obj,cons)

val = prob.solve()

fig, ax = plt.subplots()
for var, val in bluevars.items():
  print("{:10} = {:10}".format(var,val.value))
  plt.plot([val.value,val.value],[0,3])

plt.show()

结果如下所示:

如果您想“手动调整”蓝条,可以通过添加权重修改我们构建的优化问题w[i]

min   w[0]*t[0] + w[1]*t[1] + ...
s.t.  S[leftmost]  = 0
      S[i]+l[i]   <= E[i]
      E[i]        <= S+L[i]
      |(E[i]-S[i])-(L[i]+l[i])/2|_2<=t[i]

w[i]越大,所讨论的spring接近其平均长度就越重要。

根据约束最小化有序蓝条之间的平方距离

使用与上述相同的策略,假设某种已知顺序,我们可以最小化蓝色条之间的平方距离。这导致:

min   t[0] + t[1] + ...
s.t.  S[leftmost]  = 0
      S[i]+l[i]   <= E[i]
      E[i]        <= S+L[i]
      |(S[i]-S[i+1])/2|_2<=t[i]

在下面的代码中,我首先找到蓝色条的可行位置,然后假定这些映射到所需的顺序。用更准确的信息替换此启发式方法是个好主意。

#!/usr/bin/env python3

import cvxpy as cp
import networkx as nx
import matplotlib.pyplot as plt

def FindTerminalPoints(springs):
  starts = set([x[0] for x in springs.edges()])
  ends   = set([x[1] for x in springs.edges()])
  return list(starts-ends), list(ends-starts)

springs = nx.DiGraph()
springs.add_edge('a', 'b', minlen= 1, maxlen= 3)
springs.add_edge('a', 'c', minlen= 1, maxlen= 4)
springs.add_edge('a', 'f', minlen=15, maxlen=15)
springs.add_edge('b', 'c', minlen= 0, maxlen= 1)
springs.add_edge('b', 'e', minlen= 9, maxlen=11)
springs.add_edge('c', 'd', minlen= 7, maxlen=11)
springs.add_edge('d', 'e', minlen= 0, maxlen= 1)
springs.add_edge('d', 'f', minlen= 3, maxlen= 6)
springs.add_edge('e', 'f', minlen= 3, maxlen= 5)

if not nx.is_directed_acyclic_graph(springs):
  raise Exception("Springs must be a directed acyclic graph!")

starts, ends = FindTerminalPoints(springs)
if len(starts)!=1:
  raise Exception("One unique start is needed!")

if len(ends)!=1:
  raise Exception("One unique end is needed!")  

start = starts[0]
end   = ends[0]

#At this point we have what is essentially a directed acyclic graph beginning at
#`start` and ending at `end`

#Generate a variable for the position of each blue bar
bluevars = {n: cp.Variable(name=n) for n in springs.nodes()}

#Anchor the leftmost blue bar to prevent pathological solutions
cons   = [bluevars[start]==0]

#Constraint each blue bar to its range
for s,e in springs.edges():
  print("Loading edge {0}-{1}".format(s,e))
  sv   = bluevars[s]
  ev   = bluevars[e]
  edge = springs[s][e]
  cons += [sv+edge['minlen']<=ev]
  cons += [ev<=sv+edge['maxlen']]

#Find feasible locations for the blue bars. This is a heuristic for getting a
#sorted order for the bars
obj  = cp.Minimize(1)
prob = cp.Problem(obj,cons)

prob.solve()

#Now that we have a sorted order, we modify the objective to minimize the
#squared distance between the ordered bars
bar_locs = list(bluevars.values())
bar_locs.sort(key=lambda x: x.value)

dvars = [cp.Variable() for n in range(len(springs.nodes())-1)]
for i in range(len(bar_locs)-1):
  cons += [cp.norm(bar_locs[i]-bar_locs[i+1],2)<=dvars[i]]

obj  = cp.Minimize(cp.sum(dvars))
prob = cp.Problem(obj,cons)

val = prob.solve()

fig, ax = plt.subplots()
for var, val in bluevars.items():
  print("{:10} = {:10}".format(var,val.value))
  plt.plot([val.value,val.value],[0,3])

plt.show()

看起来像这样:

根据约束最小化所有蓝色条之间的平方距离

我们还可以尝试最小化蓝色条之间的所有成对平方距离。在我看来,这似乎是最好的结果。

min   t[i,j] + ...                 for all i,j
s.t.  S[leftmost]        = 0
      S[i]+l[i]         <= E[i]    for all i
      E[i]              <= S+L[i]  for all i
      |(S[i]-S[j])/2|_2 <= t[i,j]  for all i,j

看起来像这样:

#!/usr/bin/env python3

import cvxpy as cp
import networkx as nx
import matplotlib.pyplot as plt
import itertools

def FindTerminalPoints(springs):
  starts = set([x[0] for x in springs.edges()])
  ends   = set([x[1] for x in springs.edges()])
  return list(starts-ends), list(ends-starts)

springs = nx.DiGraph()
springs.add_edge('a', 'b', minlen= 1, maxlen= 3)
springs.add_edge('a', 'c', minlen= 1, maxlen= 4)
springs.add_edge('a', 'f', minlen=15, maxlen=15)
springs.add_edge('b', 'c', minlen= 0, maxlen= 1)
springs.add_edge('b', 'e', minlen= 9, maxlen=11)
springs.add_edge('c', 'd', minlen= 7, maxlen=11)
springs.add_edge('d', 'e', minlen= 0, maxlen= 1)
springs.add_edge('d', 'f', minlen= 3, maxlen= 6)
springs.add_edge('e', 'f', minlen= 3, maxlen= 5)

if not nx.is_directed_acyclic_graph(springs):
  raise Exception("Springs must be a directed acyclic graph!")

starts, ends = FindTerminalPoints(springs)
if len(starts)!=1:
  raise Exception("One unique start is needed!")

if len(ends)!=1:
  raise Exception("One unique end is needed!")  

start = starts[0]
end   = ends[0]

#At this point we have what is essentially a directed acyclic graph beginning at
#`start` and ending at `end`

#Generate a variable for the position of each blue bar
bluevars = {n: cp.Variable(name=n) for n in springs.nodes()}

#Anchor the leftmost blue bar to prevent pathological solutions
cons   = [bluevars[start]==0]

#Constraint each blue bar to its range
for s,e in springs.edges():
  print("Loading edge {0}-{1}".format(s,e))
  sv   = bluevars[s]
  ev   = bluevars[e]
  edge = springs[s][e]
  cons += [sv+edge['minlen']<=ev]
  cons += [ev<=sv+edge['maxlen']]

dist_combos = list(itertools.combinations(springs.nodes(), 2))
dvars       = {(na,nb):cp.Variable() for na,nb in dist_combos}
distcons    = []
for na,nb in dist_combos:
  distcons += [cp.norm(bluevars[na]-bluevars[nb],2)<=dvars[(na,nb)]]

cons += distcons

#Find feasible locations for the blue bars. This is a heuristic for getting a
#sorted order for the bars
obj  = cp.Minimize(cp.sum(list(dvars.values())))
prob = cp.Problem(obj,cons)

val = prob.solve()

fig, ax = plt.subplots()
for var, val in bluevars.items():
  print("{:10} = {:10}".format(var,val.value))
  plt.plot([val.value,val.value],[0,3])

plt.show()

看起来像这样: