python - 通过 L0 范数/GEKKO 中非零元素的数量约束混合整数非线性优化问题

Question

我想最小化实现一组给定的非负整数值 b 的成本，这些非负整数值 b 是从 GEKKO 中的两组非负整数变量 x、y 线性计算出来的。
如果以某种方式说明了我的问题，则 b 是 x 和 y 的约束。我的成本函数是非线性的:使用条件/最小值的二次方。除了这些标准约束之外，我还有一个约束，它要求 x 中非零元素的数量至少与 x 中的最大元素一样大(例如，L0 范数等于 LInifity 范数)。
我现在的困难是双重的，因为我对优化还很陌生，而且是 GEKKO 的新手。

我看到 GEKKO 支持 numpy 数组，这会使问题陈述变得相当简洁，但我很难让它工作 - 导致很多列表理解而不是矢量化操作。

我设法定义了 L0 范数约束，并且 GEKKO 实际使用它运行，但它未能找到解决方案。我认识到 L0 问题真的很难(例如组合)，但不知何故，解决方案很容易通过“手工”找到。我只是觉得我做错了什么。

我将不胜感激任何帮助!这是我到目前为止所做的:

from gekko import GEKKO
import numpy as np

# Setup gekko (taken from their MINLP tutorial with more iterations).
m = GEKKO()
m.options.SOLVER = 1  # APOPT is an MINLP solver
m.solver_options = ('minlp_maximum_iterations 500',
                    'minlp_max_iter_with_int_sol 10',
                    'minlp_as_nlp 0',
                    'nlp_maximum_iterations 50',
                    'minlp_branch_method 1',
                    'minlp_integer_tol 0.05',
                    'minlp_gap_tol 0.01')

# Define variables as arrays.
y = m.Array(m.Var, (7), value=1, lb=1, ub=12, integer=True)
x = m.Array(m.Var, (18), value=0, lb=0, ub=8, integer=True)

# Example of the user-defined target values b as constraints (actually input args).
m.Equation(x[2] + y[1] == 7)
m.Equation(x[12] + y[2] == 5)

# This is the L0 constraint.
# I thought using m.Array would make this a nice definition like 
#   m.Equation(np.count_nonzero(x) >= np.max(x))
# but it didn't, since these numpy functions are not supported in GEKKO.
# So I defined the following , which feels wrong:
m.Equation(m.sum([int(x_i.value > 0) for x_i in x]) >= max(x_i.value for x_i in x))

# Finally, the objective function (intermediates for readability).
k = np.array([m.min2(y_i, 3) for y_i in y])
x_cost = m.Intermediate(m.sum(x * (x + 1)))
y_cost = m.Intermediate(m.sum((k - 1) * (k + 2) + 2.5 * (y - k) * (y + k - 3)))
m.Obj(x_cost + y_cost)

# Solve.
m.solve(disp=True, debug=True)

Answer 1

Gekko 使用基于梯度的求解器，因此方程不应随迭代而改变。有一种方法可以获得与基于梯度的求解器兼容的非零变量的计数。这是另一种形式，可为您提供 x 的计数和最大值。向量。这使用 if3 , max3 , 和 min3如在 documentation on Model Building Functions with logical conditions 中找到的.

from gekko import GEKKO
import numpy as np
# Setup gekko (taken from their MINLP tutorial with more iterations).
m = GEKKO(remote=False)

# Define variables as arrays.
y = m.Array(m.Var, (7), value=1, lb=1, ub=12, integer=True)
x = m.Array(m.Var, (18), value=0, lb=0, ub=8, integer=True)

# Example of the user-defined target values b as constraints (actually input args).
m.Equation(x[2] + y[1] == 7)
m.Equation(x[12] + y[2] == 5)

# This is the L0 constraint.
#   m.Equation(np.count_nonzero(x) >= np.max(x))
eps = 0.05 # decision point for a "zero" value
count = m.sum([m.if3(x_i-eps,0,1) for x_i in x])
max_x = 0
for x_i in x:
    max_x = m.Intermediate(m.max3(max_x,x_i))
m.Equation(count >= max_x)

# Finally, the objective function (intermediates for readability).
k = np.array([m.min3(y_i, 3) for y_i in y])
x_cost = m.Intermediate(m.sum(x * (x + 1)))
y_cost = m.Intermediate(m.sum((k - 1) * (k + 2) + 2.5 * (y - k) * (y + k - 3)))
m.Minimize(x_cost + y_cost)

# Solve.
m.options.SOLVER = 3  # Initialize with IPOPT
m.solve(disp=True)

m.options.SOLVER = 1  # APOPT is an MINLP solver
m.solver_options = ('minlp_maximum_iterations 500',
                    'minlp_max_iter_with_int_sol 10',
                    'minlp_as_nlp 0',
                    'nlp_maximum_iterations 50',
                    'minlp_branch_method 1',
                    'minlp_integer_tol 0.05',
                    'minlp_gap_tol 0.01')
m.solve(disp=True)

我使用IPOPT给出初始化的非整数解，然后使用APOPT找到最佳整数解。这种方法产生了一个成功的解决方案 x

>>> x
array([[0.0], [0.0], [0.0], [0.0], [0.0], [0.0], [0.0], [0.0], [0.0],
       [0.0], [0.0], [0.0], [0.0], [0.0], [0.0], [0.0], [0.0], [0.0]],
      dtype=object)

和

>>> y
array([[1.0], [7.0], [5.0], [1.0], [1.0], [1.0], [1.0]], dtype=object)

您可能不打算为 x 提供零解决方案因此您可能需要添加约束，例如 m.Equation(count>=n)或更改 eps=0.05找到一个非零值或在局部最小值处将求解器推离零。与 m.Equation(count>=n)和 eps=0.5它找到了更好的解决方案:

n=3, obj=63
n=5, obj=52

这是 x的解决方案和 y当 obj=52 (找到最佳解决方案)。

>>> x
array([[0.0], [1.0], [4.0], [0.0], [0.0], [0.0], [0.0], [0.0], [0.0],
       [0.0], [1.0], [0.0], [2.0], [0.0], [0.0], [0.0], [1.0], [0.0]],
      dtype=object)

>>> y
array([[1.0], [3.0], [3.0], [1.0], [1.0], [1.0], [1.0]], dtype=object)

您可能需要调整 eps用于 if3当一个值被计为非零时进行调整或调整 minlp_integer_tol 0.05这是确定整数容差的求解器选项。这是带有附加不等式约束的最终脚本，可提供最佳解决方案。

from gekko import GEKKO
import numpy as np
# Setup gekko (taken from their MINLP tutorial with more iterations).
m = GEKKO(remote=False)

# Define variables as arrays.
y = m.Array(m.Var, (7), value=1, lb=1, ub=12, integer=True)
x = m.Array(m.Var, (18), value=0, lb=0, ub=8, integer=True)

# Example of the user-defined target values b as constraints (actually input args).
m.Equation(x[2] + y[1] == 7)
m.Equation(x[12] + y[2] == 5)

# This is the L0 constraint.
#   m.Equation(np.count_nonzero(x) >= np.max(x))
eps = 0.5 # threshold for a "zero" value
count = m.sum([m.if3(x_i-eps,0,1) for x_i in x])
max_x = 0
for x_i in x:
    max_x = m.Intermediate(m.max3(max_x,x_i))
m.Equation(count >= max_x)

# Finally, the objective function (intermediates for readability).
k = np.array([m.min3(y_i, 3) for y_i in y])
x_cost = m.Intermediate(m.sum(x * (x + 1)))
y_cost = m.Intermediate(m.sum((k - 1) * (k + 2) + 2.5 * (y - k) * (y + k - 3)))
m.Minimize(x_cost + y_cost)

m.Equation(count>=5)

# Solve.
m.options.SOLVER = 3  # Initialize with IPOPT
m.solve(disp=True)

m.options.SOLVER = 1  # APOPT is an MINLP solver
m.solver_options = ('minlp_maximum_iterations 500',
                    'minlp_max_iter_with_int_sol 10',
                    'minlp_as_nlp 0',
                    'nlp_maximum_iterations 50',
                    'minlp_branch_method 1',
                    'minlp_integer_tol 0.05',
                    'minlp_gap_tol 0.01')
m.solve(disp=True)

您也许可以调整 n或一些求解器选项以获得更好的解决方案。我希望这有助于为您指明正确的方向。