python - 如何解决 Python Gekko 中的积分问题？

Question

如何在 Python Gekko 的优化模型中包含变量 x 的积​​分？我发现 m.vsum(x) 是时间范围内的求和，但不考虑可变时间步长。变量 x.dt() 的导数是 Gekko 内置的，但是积分呢？

from gekko import GEKKO
import numpy as np
m = GEKKO()
m.time = np.linspace(0,2,5)
m.options.IMODE=4; m.options.NODES=3
x = m.Var(5)
m.Equation(x.dt()==-x)
#y = m.Var(0)
#m.Equation(y=m.integral(x))
m.options.SOLVER=1
m.solve()
print(x.value)

结果

[5.0, 3.0434782609, 1.8525519849, 1.1276403386, 0.68638977133]

Answer 1

Gekko 求解微分方程和代数方程。您可以通过定义一个导数等于 x 的新变量 y 来包含积分。

from gekko import GEKKO
import numpy as np
m = GEKKO()
m.time = np.linspace(0,2,5)
m.options.IMODE=4; m.options.NODES=6
x = m.Var(5)
m.Equation(x.dt()==-x)
y = m.Var(0)
m.Equation(y.dt()==x)
m.solve()
print(np.round(y.value,5))

结果是 [0., 1.96735, 3.1606, 3.88435, 4.32332] 即 x 在时间点 [0, 0.5, 1.0, 1.5] 的积分, 2.0]。 y 的初始条件为零，但它可以是非零的，例如比例积分微分 (PID) Controller ，其中积分项可能在每个循环中累加。

在 Gekko v0.2.6 (pip install gekko --upgrade) 中，将有一个新的 integral 函数(感谢您的提问)。您的注释行是正确的，但不要忘记方程表达式中的双等号。

m.Equation(y==m.integral(x))

这是数值解和精确解的比较。

from gekko import GEKKO
import numpy as np
m = GEKKO()
m.time = np.linspace(0,2,5)
m.options.IMODE=4; m.options.NODES=6
x = m.Var(5)
m.Equation(x.dt()==-x)
y = m.Var(0)
m.Equation(y==m.integral(x))
z = m.Var(0)
m.Equation(z.dt()==x)
m.options.SOLVER=1
m.solve()

# results
print('Numerical x')
print(np.round(x.value,6))
print('Exact x')
print(np.round(5*np.exp(-m.time),6))
print('Numerical Integral y')
print(np.round(y.value,6))
print('Numerical Integral z End Point')
print(np.round(z.value[-1],6))
print('Exact Integral')
print(np.round(-5*(np.exp(-2)-np.exp(0)),6))

结果有 1e-6 个差异，因为它们是仅在请求的时间点计算的数值解。附加时间点的解决方案更准确，但解决方案速度较慢(尤其是对于大型问题)。

Numerical x
[5.       3.032654 1.839398 1.115651 0.676677]
Exact x
[5.       3.032653 1.839397 1.115651 0.676676]
Numerical Integral y
[0.       1.967346 3.160602 3.884349 4.323323]
Numerical Integral z End Point
4.323323
Exact Integral
4.323324