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

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

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

我现在的困难是双重的，因为我对优化还很陌生，而且是 GEKKO 的新手。

1. 我看到 GEKKO 支持 numpy 数组，这会使问题陈述变得相当简洁，但我很难让它工作 - 导致大量列表推导而不是矢量化操作。
2. 我设法定义了 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)

解决方法

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

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,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]],dtype=object)

和 y

>>> y
array([[1.0],[7.0],[5.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],[4.0],[2.0],dtype=object)

>>> y
array([[1.0],[3.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,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,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_gap_tol 0.01')
m.solve(disp=True)

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