如何在 GEKKO 中模拟不随时间而随体积变化的微分方程？

如何解决如何在 GEKKO 中模拟不随时间而随体积变化的微分方程？

我知道正常的方法。例如...

m.Equation(x.dt() == v)

但是对于像 PFR 这样的东西呢？假定处于稳态，但流速取决于反应器的体积（或长度）。如果我还想优化模型中的某些内容，我觉得 IMODE = 3 是最好的。我认为它看起来像这样。

m.Equation(Fa.dVr() == -ra)

有什么帮助吗？我只是将 dt() 视为 dVr() 吗？并将 m.time 声明为 Volume 数组？它是否被视为依赖于体积的动态系统？

解决方法

我刚刚测试了这个。下面的结果与我在 odeint 中执行时产生的结果相同

from gekko import GEKKO
import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import odeint



k = .00384
Keq = 4
KB = .46
KW = 3.20
Ca_0 = 1.5
Cb_0 = 1.5
Cc_0 = 0
Cd_0 = 0
Fa_0 = 1 # mol/h
v = Fa_0/Ca_0
Fb_0 = Cb_0*v
Fc_0 = 0
Fd_0 = 0

m = GEKKO()
m.time = np.linspace(0,4000,1000)

Fa = m.Var(Fa_0)
Fb = m.Var(Fb_0)
Fc = m.Var(Fc_0)
Fd = m.Var(Fd_0)

Ca = m.Intermediate(Fa/v)
Cb = m.Intermediate(Fb/v)
Cc = m.Intermediate(Fc/v)
Cd = m.Intermediate(Fd/v)

r = m.Intermediate(k*(Ca*Cb - Cc*Cd/Keq)/\
                  ( 1 + KB*Cb + KW*Cd))

m.Equation(Fa.dt() == -r)
m.Equation(Fb.dt() == -r)
m.Equation(Fc.dt() ==  r)
m.Equation(Fd.dt() ==  r)

m.options.IMODE = 4
m.solve(disp=False)

plt.plot(m.time,Fa.value,label='Fa')
plt.plot(m.time,Fb.value,label='Fb')
plt.plot(m.time,Fc.value,label='Fc')
plt.plot(m.time,Fd.value,label='Fd')

from gekko import GEKKO
import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import odeint



k = .00384
Keq = 4
KB = .46
KW = 3.20
Ca_0 = 1.5
Cb_0 = 1.5
Cc_0 = 0
Cd_0 = 0
Fa_0 = 1 # mol/h
v = Fa_0/Ca_0
Fb_0 = Cb_0*v
Fc_0 = 0
Fd_0 = 0

def f(F,W):
    Fa,Fb,Fc,Fd = F
    Ca,Cb,Cc,Cd = Fa/v,Fb/v,Fc/v,Fd/v
    r = k*(Ca*Cb - Cc*Cd/Keq)/( 1 + KB*Cb + KW*Cd)
    
    return -r,-r,r,r
W = np.linspace(0,1000)
Fa1,Fb1,Fc1,Fd1 = odeint(f,[Fa_0,Fb_0,Fc_0,Fd_0],W).T

plt.plot(W,Fa1,label='Fa')
plt.plot(W,label='Fb')
plt.plot(W,label='Fc')
plt.plot(W,Fd1,label='Fd')