问题描述
最近,我用GEKKO建立了一个小模型.它包含一个实际上随时间变化的参数.我该如何实施?我尝试使用if3
,但它给出了一个错误.
Recently, I have built a small model with GEKKO. It contains a Parameter which actually changes with time. How can I implement that? I tried using if3
, but it gives an error.
这是MWE:
#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Started on 10-08-2019
@author: rotton
"""
import numpy as np
import matplotlib.pyplot as plt
from gekko import GEKKO
#Initialize Model
m = GEKKO(remote=False)
# Parameters
k_1 = m.Param(value = 0.19)
f_1 = m.Param(value = 29.0)
V_liq = m.Param(value = 159.0)
q_in = m.Param(value = 2.5)
X_in = m.Param(value = 271.77)
Y_in = m.Param(value = 164.34)
X = m.Var(value = 11.55)
Y = m.Var(value = 11.55*0.2)
rho_1 = m.Intermediate(k_1*X)
q_prod = m.Intermediate(0.52*f_1*X)
m.time = np.arange(0,5,1/12)
m.Equations([X.dt() == q_in/V_liq*(X_in - X) - rho_1, \
Y.dt() == q_in/V_liq*(Y_in - Y)])
#Dynamic simulation
m.options.IMODE = 4
m.solve(disp=False)
plt.plot(m.time, X.value)
plt.xlabel('time')
plt.ylabel('X')
plt.show()
我尝试了以下操作:
q_in = m.if3(m.time - 2, 0, 2.5)
,因此q_in
最初将为0,在time = 2
处变为2.5.但是我收到以下错误:
so that q_in
would be 0 initially, and become 2.5 at time = 2
.But I get the following error:
File "/usr/local/lib/python3.7/site-packages/gekko/gekko.py", line 1838, in solve
raise Exception(apm_error)
Exception: @error: Equation Definition
Equation without an equality (=) or inequality (>,<)
(((1-int_v5))*([-2.-1.91666667-1.83333333-1.75-1.66666667-1.58333333
STOPPING...
您是否知道我该如何实现?实际上,此变量在0到60之间跳跃了几次,而我在CSV文件中有可用的时间点.理想情况下,我可以创建一个循环,该循环将在每次迭代时检查q_in
是否需要更改,并相应地覆盖当前值.
Do you have an idea how I can achieve this? Actually, this variable jumps several times between 0 and 60, and I have the time points available in a CSV file. Ideally, I could create a loop that would check at every iteration if it's time for q_in
to change, and overwrite the current value accordingly.
推荐答案
您可以从CSV中读取输入,并将时变值分配给q_in.value
,无论是在参数初始化期间(请参见示例1)还是在循环,每次积分间隔值都会改变(请参见示例2).示例1和2都产生以下结果,但是示例1更快.
You can read the input from a CSV and assign the time-varying values to q_in.value
either during Parameter initialization (see Example #1) or else in a loop where the value changes each time integration interval (see Example #2). Examples 1 and 2 both produce the following result but Example 1 is faster.
如果时间范围很长,使用选项m.options.IMODE=7
的示例1可能也会更快. IMODE=7
使用顺序求解方法而不是同时求解方法.
Example 1 may also be faster with option m.options.IMODE=7
if you have a very long time horizon. IMODE=7
uses a sequential solution method instead of a simultaneous solution method.
示例1
import numpy as np
import matplotlib.pyplot as plt
from gekko import GEKKO
t = np.arange(0,5,1/12)
step = [0 if z<2 else 2.5 for z in t]
m = GEKKO(remote=False)
k_1 = m.Param(value = 0.19)
f_1 = m.Param(value = 29.0)
V_liq = m.Param(value = 159.0)
q_in = m.Param(value = step)
X_in = m.Param(value = 271.77)
Y_in = m.Param(value = 164.34)
X = m.Var(value = 11.55)
Y = m.Var(value = 11.55*0.2)
rho_1 = m.Intermediate(k_1*X)
q_prod = m.Intermediate(0.52*f_1*X)
m.time = t
m.Equations([X.dt() == q_in/V_liq*(X_in - X) - rho_1, \
Y.dt() == q_in/V_liq*(Y_in - Y)])
m.options.IMODE = 4
m.solve(disp=False)
plt.plot(m.time,q_in.value,label=r'$q_{in}$')
plt.plot(m.time, X.value,label='X')
plt.plot(m.time, Y.value,label='Y')
plt.legend()
plt.xlabel('time')
plt.show()
示例2
import numpy as np
import matplotlib.pyplot as plt
from gekko import GEKKO
t = np.arange(0,5,1/12)
m = GEKKO(remote=False)
k_1 = m.Param(value = 0.19)
f_1 = m.Param(value = 29.0)
V_liq = m.Param(value = 159.0)
q_in = m.Param()
X_in = m.Param(value = 271.77)
Y_in = m.Param(value = 164.34)
X = m.Var(value = 11.55)
Y = m.Var(value = 11.55*0.2)
rho_1 = m.Intermediate(k_1*X)
q_prod = m.Intermediate(0.52*f_1*X)
m.time = [t[0],t[1]]
m.Equations([X.dt() == q_in/V_liq*(X_in - X) - rho_1, \
Y.dt() == q_in/V_liq*(Y_in - Y)])
m.options.IMODE = 4
# store Xs and Ys for plotting
for i in range (1,len(t)):
q_in.value = 0 if t[i]<2 else 2.5
m.solve(disp=False)
if i==1:
Xs = [X.value[0]]
Ys = [Y.value[0]]
Xs.append(X.value[1])
Ys.append(Y.value[1])
step = [0 if z<2 else 2.5 for z in t]
plt.plot(t,step,label=r'$q_{in}$')
plt.plot(t, Xs,label='X')
plt.plot(t, Ys,label='Y')
plt.legend()
plt.xlabel('time')
plt.show()
如果需要使q_in
依赖于某些变量的值,则可以使用m.if3
函数.但是,这是一个更具挑战性的问题,因为m.if3
函数将问题转换为混合整数非线性编程形式,可能需要更长的时间才能解决.这是q_in=0
时X>8
和q_in=2.5
时X<=8
的示例.但是,它对我来说并不收敛.我不知道为什么,我需要做一些额外的挖掘工作,但是我希望您能拥有它,以防它对您有用.
If you need to make q_in
dependent on the value of some of your variables then you can use the m.if3
function. However, this is a more challenging problem to solve because the m.if3
function converts the problem into a Mixed Integer Nonlinear Programming form that may take longer to solve. Here is an example where q_in=0
when X>8
and q_in=2.5
when X<=8
. However, it didn't converge for me. I'm not sure why and I'd need to do some additional digging but I though you'd like to have it in case it works for you.
import numpy as np
import matplotlib.pyplot as plt
from gekko import GEKKO
m = GEKKO(remote=False)
k_1 = m.Param(value = 0.19)
f_1 = m.Param(value = 29.0)
V_liq = m.Param(value = 159.0)
X_in = m.Param(value = 271.77)
Y_in = m.Param(value = 164.34)
X = m.Var(value = 11.55,name='X')
Y = m.Var(value = 11.55*0.2,name='Y')
rho_1 = m.Intermediate(k_1*X)
q_prod = m.Intermediate(0.52*f_1*X)
q_in = m.if3(8-X, 0.0, 2.5)
m.time = np.arange(0,5,1/12)
m.Equations([X.dt() == q_in/V_liq*(X_in - X) - rho_1, \
Y.dt() == q_in/V_liq*(Y_in - Y)])
m.options.IMODE = 6
m.options.SOLVER = 1
m.solve(disp=True)
plt.plot(m.time,q_in.value,label=r'$q_{in}$')
plt.plot(m.time, X.value,label='X')
plt.plot(m.time, Y.value,label='Y')
plt.legend()
plt.xlabel('time')
plt.show()
还有一个这里没有其他示例可以解决随着时间的ODE问题使用Gekko进行各种输入.
There are also a few other examples here on solving ODEs with time-varying inputs with Gekko.
这篇关于如何定义时间相关的离散参数?的文章就介绍到这了,希望我们推荐的答案对大家有所帮助,也希望大家多多支持!