使用 Gekko 和 Python 拟合数据的数值 ODE 解决方案

Question

使用 Gekko 为数据拟合数值 ODE 解。

大家好！ 我想知道，是否可以使用 GEKKO 拟合 ODE 的系数。 我尝试复制 example given here 失败。

这就是我想出的（但有缺陷——我也许应该提到我的数学技能很差）：

import numpy as np
from gekko import GEKKO

tspan = [0, 0.1, 0.2, 0.4, 0.8, 1]
Ca_data = [2.0081,  1.5512,  1.1903,  0.7160,  0.2562,  0.1495]

m = GEKKO(remote=False)

t = m.Param(value=tspan)
m.time = t
Ca_m = m.Param(value=Ca_data)

Ca = m.Var()

k = m.FV(value=1.3)
k.STATUS = 1

m.Equation( Ca.dt() == -k * Ca)

m.Obj( ((Ca-Ca_m)**2)/Ca_m )

m.options.IMODE = 2
m.solve(disp=True)
print(k.value[0]) #2.58893455 is the solution

有人可以帮我吗？ 非常感谢你， 马丁

（这是我在这里的第一篇文章——如果我做了不恰当的事情，请保持温和。）

Answer 1

您的解决方案很接近，但您需要：

• 更多节点（默认=2）以提高准确性。 Gekko 只添加您定义的点。见additional information on collocation。
• 将Ca 定义为m.CV() 以使用内置错误模型而不是m.Var() 和m.Obj with NODES>=3。否则，每个搭配区间的内部节点也与测量值匹配，这给出了一个稍微错误的答案。
• 设置EV_TYPE=2 使用平方误差。绝对值目标 EV_TYPE=1（默认）给出正确但略有不同的答案。

import numpy as np
from gekko import GEKKO

m = GEKKO(remote=False)
m.time = [0, 0.1, 0.2, 0.4, 0.8, 1]
Ca_data = [2.0081,  1.5512,  1.1903, 0.7160,  0.2562,  0.1495]
Ca = m.CV(value=Ca_data); Ca.FSTATUS = 1 # fit to measurement
k = m.FV(value=1.3); k.STATUS = 1        # adjustable parameter
m.Equation(Ca.dt()== -k * Ca)            # differential equation

m.options.IMODE = 5   # dynamic estimation
m.options.NODES = 5   # collocation nodes
m.options.EV_TYPE = 2 # squared error
m.solve(disp=True)    # display solver output
print(k.value[0])     # 2.58893455 is the curve_fit solution

解决方案是k=2.5889717102。绘图显示与测量值的匹配。

import matplotlib.pyplot as plt  # plot solution
plt.plot(m.time,Ca_data,'ro')
plt.plot(m.time,Ca.value,'bx')
plt.show()

有additional tutorials和course material关于微分和代数方程模型的参数估计。