后续问题:矩阵形式的GEKKO优化

Question

这是我之前发布的问题的后续问题：

矩阵形式的GEKKO优化

我需要再添加一个跟踪"库存"("Inv"(的约束，它跟踪总和(q[I，：]-q[：，I](。"Inv"将是一个4X1列矢量。我尝试了以下方法：

m = GEKKO(remote=False)
q = m.Array(m.Var,(4,4),lb=0,ub=10)
for i in range(4):
for j in range(4):
if j<=i:
q[i,j].upper=0 # set upper bound = 0
def profit(q):
profit = np.sum(q.flatten() * pmx.flatten())
return profit
Inv[0]=0
for i in range(4):
m.Equation(np.sum(q[i,:])<=10)
m.Equation(np.sum(q[:,i])<=8)
m.Equation(I[i] = I[i-1] + (np.sum(q[i,:]) - np.sum(q[:,i]))) # New Line 1a inventory 
Inv[i] = Inv[i-1] + (np.sum(q[i,:]) - np.sum(q[:,i])) # New Line 1b inventory. Keep either 1a or 1b 
m.Equation(Inv[i] <= 15) # New Line 2 inventory constraint
m.Equation(Inv[4] = 0) # New Line 3 ending inventory should be equal to starting inventory
m.Maximize(profit(q))
m.solve()
print(q)
qr = np.array([[q[i,j].value[0] for j in range(4)] for i in range(4)])
Ir = np.array([Inv[i].value[0] for i in range(4)]) #New Line 4

错误：

1a。添加新行1a："关键字不能是表达式">

1b。用1b替换新线1a：没有问题(但是，我不确定GEKKO是否会跟踪I。此外，我需要定义"I"，即"q"的方式……不确定如何定义(。替换=注释掉1a，然后用1b运行代码。

3. 新行2：错误="int类型的对象没有len(("；但类型(I(显示为ndarray。(保留了新的1a和1b号线，然后增加了新的2号线(

4. 新建第3行：错误="关键字不能是表达式"(保留第32行，然后添加第3行(

5. 第4行：错误"'numy.ndarray'对象没有属性'value'"[删除了第3行和第4行。这是有道理的，因为如果我不能在模型中捕捉"库存"，那么它就没有价值属性(

问题：1.我对存货的定义是否正确？

2. 如果是，可以在当前型号规范中进行，还是需要完全不同的配方？如果是的话，你能介绍一下会是什么吗？

3. 在GEKKO网站上发布的各种视频中，有没有一个特定的视频我应该查看以获取更多信息？我在想DO视频，但我不认为这是一个动态优化问题(因为我没有试图优化最佳路径(。。。。

再次感谢您的帮助，

----更新5/10也尝试过：

Inv = m.SV()
for i in range(4):
m.Equation(np.sum(q[i,:])<=10)
m.Equation(np.sum(q[:,i])<=8)
#m.Equation(I[i] = I[i-1] + (np.sum(q[i,:]) - np.sum(q[:,i])))
m.Equation(m.Inv.dt() == m.Inv + (np.sum(q[i,:]) - np.sum(q[:,i])))
#I[i] = I[i-1] + (np.sum(q[i,:]) - np.sum(q[:,i])
m.Equation(m.Inv <= 15)
#m.Equation(I[4] = 0)
m.Maximize(profit(q))

新错误："GEKKO"对象没有属性"Inv">

Answer 1

实现这一点的一种方法是从Inv[0]=0的零库存开始，然后使用gekko变量Inv[1:4]跟踪库存金额。关于构建模型的几个提示：

• 对相等约束使用双等号
• 你可以定义一个临时变量，比如I=Inv[2]+Inv[3]，但它不会是gekko变量
• 您可能还想查看那些显式计算的中间变量。这样可以加快计算速度
• 我在动态优化课程中推荐本教程

import numpy as np
import scipy.optimize as opt
from gekko import GEKKO
p= np.array([4, 5, 6.65, 12]) #p = prices
pmx = np.triu(p - p[:, np.newaxis]) #pmx = price matrix, upper triangular
m = GEKKO(remote=False)
q = m.Array(m.Var,(4,4),lb=0,ub=10)
# only upper triangular can change
for i in range(4):
for j in range(4):
if j<=i:
q[i,j].upper=0 # set upper bound = 0
def profit(q):
profit = np.sum(q.flatten() * pmx.flatten())
return profit
Inv = m.Array(m.Var,5,lb=0,ub=15)
Inv[0].upper = 0 # start with 0 inventory
for i in range(4):
m.Equation(np.sum(q[i,:])<=10)
m.Equation(np.sum(q[:,i])<=8)
# track inventory    
m.Equation(Inv[i+1]==Inv[i] + (m.sum(q[i,:])-m.sum(q[:,i]))) 
m.Equation(Inv[4] == 0) # Use double == sign, not needed in loop
m.Maximize(profit(q))
m.solve()
print(q)
# convert to matrix form
qr = np.array([[q[i,j].value[0] for j in range(4)] for i in range(4)])
for i in range(4):
rs = qr[i,:].sum()
print('Row sum ' + str(i) + ' = ' + str(rs))
cs = qr[:,i].sum()
print('Col sum ' + str(i) + ' = ' + str(cs))    
Ir = np.array([Inv[i].value[0] for i in range(4)])
print(Ir)