scipy.optimize 具有多个边界、约束和连续字段答案

【问题标题】：scipy.optimize with multiple bounds, constraints and continuous fieldsscipy.optimize 具有多个边界、约束和连续字段
【发布时间】：2018-10-12 14:40:10
【问题描述】：

我想优化CHP plant 在请求的电源配置文件上的操作。因此，我定义了一个功率配置文件，热电联产工厂应尽可能遵循该配置文件。必须应用多个界限和约束来表示 CHP 工厂的实际运行。例如，这包括 CHP 可以打开或关闭，并且在打开时，其功率调制只能设置为特定的百分比范围。

下面是一个简短的解释的最小工作示例：

import scipy.optimize as opt
import numpy as np

x = np.arange(200)  # dummy x vector
poly_profile = np.array(  # 7th degree polynome fit of profile
    [-2.14104340e-11,  1.85108903e-08, -6.66697810e-06,  1.29239710e-03,
     -1.45110876e-01,  9.40324129e+00, -3.24548750e+02,  4.60006330e+03])
poly_fun = np.poly1d(poly_profile)  # make poly fun
profile = poly_fun(x[65:196])
x0 = np.zeros_like(profile)  # all zeros as starting values

def optifun(x, profile):  # define minimization fun
    return - np.sum(profile * x)

bnds_hi = opt.Bounds(0.3, 1)  # upper bounds
bnds_lo = opt.Bounds(0, 0)  # lower bounds

res = opt.minimize(
    optifun, x0, args=(profile), bounds=bnds_hi,
    constraints={'type': 'eq', 'fun': lambda x:  np.sum(x*40) - 2000},
    method='SLSQP')
plt.plot(res.x)
plt.plot(profile)

所以这些是我想要使用的界限：

(x == 0) or (0.3 <= x <= 1)，对于数组中的任意值 x
这意味着总 CHP 功率的调制度 x 可以是 0（关闭）或>0.3 和<= 1。但我可以指定下限 OR 上限。仅指定上限就不可能“关闭 CHP”，而将下限设置为 bnds_lo = opt.Bounds(0, 1)
将使热电联产工厂能够在不现实的运行点（功率调制的 0% 到 30% 之间）运行。
有没有办法让这个工作在最小工作示例中指定的范围内工作？具体来说：我可以同时设置两种边界，比如bounds=[bnds_lo, bnds_hi]吗？
我想这是一个混合整数线性规划问题，但 COBYLA 或 SLSQP 不应该能够处理这个问题吗？如果没有：是否有任何解决方法？

以及我想使用的约束：

np.sum(x*40) - 450
将热输出限制为某些热存储容量。这里 40 是热输出功率，而 450 是剩余存储容量。这很容易实现。
限制 CHP 工厂的启动次数。作为一个例子，我们假设
```
bnds_lo = opt.Bounds(0, 1)  # lower bounds
res = opt.minimize(
    optifun, x0, args=(profile), bounds=bnds_lo,
    constraints={'type': 'eq', 'fun': lambda x:  np.sum(x*40) - 1000},
    method='SLSQP')
```
这导致热电联产工厂运行 3 个周期。有没有办法限制这个？我正在考虑添加一个特定的约束函数，该函数在前导 0 之后计算正差异，但我无法完成类似的工作（例如，由于大多数 x 不完全为 0，因为边界设置为 @987654332 @.但是其他问题也可能是原因）...
设置 CHP 工厂的最短连续运行时间。这意味着至少有 5 个连续的x != 0 应该是有利的。我想过尝试与上一点类似的东西（限制启动次数），但也无法找到有用的东西。这是迄今为止最不重要的问题。

为了解决这些问题，我也尝试使用 scipy.optimize.LinearConstraings 和 NonlinearConstraings 但是method='trust-constr' 需要一个 jac（据我在 github 上阅读，这似乎是一个错误），因此我无法使其工作。

有什么方法可以让我完成这项工作吗？尤其是指定多个边界很重要。

提前致谢！

真诚地，斯科蒂

【问题讨论】：

由于问题的离散性（开关决策），这些类型的电力调度模型通常使用混合整数规划模型来解决。边界(x = 0) or (0.3 <= x <= 1) 通常被描述为 x 是一个半连续变量。
感谢您的回答。我想避免使用像 pyomo 这样的 MILP 工具，以减少我的程序对其他模块的依赖。在一个更大的仿真程序中，每隔几步就会调用这种优化，并将用于设置仿真中的一些变量。因此，我真的很想避免使用过于复杂的优化来保持良好的性能。如果没有直接的方法来实现离散边界，我仍然可以在优化后手动裁剪它们。这当然是最坏的情况，但还是可以的，因为每隔几步就会执行一次优化，因此会出错
会很小。但我也无法弄清楚如何限制启动次数。关于如何做到这一点的任何想法？
只关注约束(x == 0) or (0.3 <= x <= 1)：scipy 最小化器无法处理这种类型的约束。要为此使用 scipy，您必须进行两次最小化，一次使用 x == 0，再次使用 0.3 <= x <= 1，然后选择最佳答案。
好的，所以我想我只需要使用 0 <= x <= 1 或 milp 优化器来完成。谢谢！对处理约束有什么建议吗？

标签： python numpy optimization scipy

【解决方案1】：

profile * x0 在您的代码中给出
“ValueError：操作数无法与形状 (131,) (200,) 一起广播”。

只是猜测，x_t 是产品onoff_t * xon_t
onoff_t = 0 或 1
和0.3 <= xon_t <= 1 在t 中的每个0 .. T ?
IE。对于T = 5，有 2^5 个可能的onoff 序列，00000 00001 00010 .. 11111 ?

如果是这样，最大化sum 0:T w_t * onoff_t * xon_t 具有固定权重函数w_t 是微不足道的：
其中w_t <= 0：onoff_t = 0，关闭
其中w_t > 0：onoff_t = 1，开，xon_t = 1，最大值。
所以这不是你的问题——请澄清。

如果onoff_t 被进一步限制为只能切换两次，0...1...0...，那么可能的序列数量足够小，可以全部尝试，顺其自然：

def pulse_generator( T=200, minwidth=5 ):
    """ -> arrays of T floats, 0... 1... 0... """
    for t0 in xrange( 1, T ):
        for t1 in xrange( t0 + minwidth, T ):
            pulse = np.zeros( T )
            pulse[t0:t1] = 1
            yield pulse

for pulse in pulse_generator( T ):
    print "pulse:", pulse
    optimize myfunction( pulse * xon ), 0.3 <= xon <= 1

切换4次，0... 1... 0... 1... 0...，类似。（对于给定的T，有多少这样的脉冲？见维基百科Stars and bars—— 太棒了。）

补充：我不是专家，但也不是开关 bang-bang control 非常对微小的变化敏感，早一点还是晚一点？一个程序（mer）可能会花费大量时间在噪音中进行抖动。两个阶段怎么样，粗网格然后细网格 -

将时间 0:T 分成 10 段，运行所有 2^10 = 1024 个开关序列
1a。仔细看看最好的——有什么规律吗？
将它们的边缘移动半步，T / 20。

另请参阅：谷歌“离散优化”多重网格......和 Grid search.

【讨论】：

嗨，丹尼斯。非常感谢你的帮助。是的，既然您给了我一个将x 分隔为onoff_t 和xon_t 的示例，这似乎微不足道。但是当我试图自己解决它时，我并没有考虑到这种可能性。 :) 脉冲发生器也很有帮助。但是我不能让它适用于任意数量的脉冲，其中只有最大脉冲数是有限的，但允许使用更少的脉冲。例如最多 4 个脉冲，但最佳解决方案只使用 2 个脉冲。我应该就这个主题开一个新问题吗？
profile * x0 不会在我的代码中引发任何错误。好像您使用的是 python 2.x。在 Python 3 和 2.x 之间使用 numpy 数组可能存在一些差异。因为profile = poly_fun(x[65:196])x0 = np.zeros_like(profile) 行应该清楚地避免profile 和x0 的任何形状不匹配。
@Scotty1- 假设您如上所述定义了pulsegen2 和pulsegen3。运行它们的一种方法是：for pulse in list(pulesegen2) + list(pulsegen3): ... 一种为大型生成器使用更少内存的更好方法是itertools.chain，请参阅Functional Programming Howto。任意数量的脉冲都是代码高尔夫，首先执行 2 3 4。此外，“固定权重函数......微不足道”？
感谢您提供更多信息。是的，bang-bang control 并不是处理控件的好方法，但在许多领域它仍然是一种常见的控制策略。此外，正如您所建议的那样，将“bang-bang”脉冲发生器与连续变量叠加听起来像是一种可接受的解决方法，因为这种情况下的优化问题相当小，并且结果被输入到一个主要由非刚性 ODE 组成的系统中惯性大，导致整个系统对优化结果突变的敏感性低。
关于函数式编程：在我的理解中，函数式编程意味着编写一个通用的脉冲发生器函数，该函数可以在脉冲发生器中重复用于任意数量的脉冲，而无需编写特定数量的脉冲发生器的脉冲。 Itertools.chain 似乎仍然是处理特定数量脉冲的脉冲发生器的好方法。