使用 Gekko 进行整数优化的问题答案

【问题标题】：Issue with integer optimization with Gekko使用 Gekko 进行整数优化的问题
【发布时间】：2021-01-02 18:54:22
【问题描述】：

我尝试了这段代码，结果是(2,2)。显然这不是真的，因为(0,4) 是最好的解决方案。你能告诉我问题出在哪里吗？

import numpy as np
from gekko import GEKKO
m = GEKKO()

x=[m.Var(lb=0,ub=4,integer=True) for i in range(2)]

m.options.SOLVER = 1

#Equations

m.Equation(m.sum(x)>=4)

#Objective

m.Obj(x[0]*x[1])
m.solve(disp=False)
print(x[0].value)
print(x[1].value)

【问题讨论】：

它使用 m.Equation(x[0]*x[1]<=1e-5) 而不是 m.Obj(x[0]*x[1])。
@JunhoPark 为什么要添加约束？它不在模型中。由求解器来最小化目标函数

标签： python numpy optimization integer gekko

【解决方案1】：

所有求解器（m.options.SOLVER 和 1=APOPT、2=BPOPT、3=IPOPT）通过初始猜测 (0,0) 和其他初始猜测找到局部最大值 (2,2)。在(0,4) 和(4,0) 有两个局部最优值。这是一个非凸优化问题。求解器是非线性凸优化求解器，因此它们“应该”找到局部解之一。如果您旋转图形以与x0=x1 线对齐，这个问题让我想起了saddle point problem。如果不检查不确定的 Hessian，优化器可能会遇到鞍点问题。

更好的初始猜测有助于求解器找到最优解，例如x0=1 和x1=3。此外，变量之一的上限<1.999 或>2.001 的下限也可以找到最佳值。约束限制了搜索区域，使其成为凸优化问题并避免了鞍点。

from gekko import GEKKO
m = GEKKO()

x = m.Array(m.Var,2,value=1,lb=0,ub=4,integer=True)
x[0].value = 1; x[1].value = 3 # initial guesses
#x[0].UPPER=1
m.Equation(m.sum(x)>=4)
m.Minimize(x[0]*x[1])

m.options.SOLVER = 1
m.solve(disp=True)
print(x[0].value)
print(x[1].value)

我创建了优化问题的等高线和 3D 图。由于问题是非凸的，您可能希望切换到多启动方法或使用非凸优化求解器。

import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
from matplotlib import cm

# Design variables at mesh points
x0 = np.arange(0.0, 4.01, 0.1)
x1 = np.arange(0.0, 4.01, 0.1)
x0,x1 = np.meshgrid(x0,x1)

# Equation sum(x)>4
xsum = x0+x1

# Objective Minimize(x0*x1)
xobj = x0*x1

# Create a contour plot
plt.figure()
CS = plt.contour(x0, x1, xobj)
plt.clabel(CS, inline=1, fontsize=10)
CS = plt.contour(x0, x1, xsum,[3.6,3.8,4.0],\
                 colors='k',linewidths=[0.5,1,4])
plt.clabel(CS, inline=1, fontsize=10)
plt.xlabel('x0')
plt.ylabel('x1')
plt.savefig('contour.png')

# Create 3D plot
fig = plt.figure()
ax = fig.gca(projection='3d')
surf = ax.plot_surface(x0, x1, xobj, \
            rstride=1, cstride=1, cmap=cm.coolwarm, \
            vmin = 0, vmax = 10, linewidth=0, \
            antialiased=False)
plt.show()

【讨论】：