如何定义时间相关的离散参数

如何定义时间相关的离散参数

本文介绍了如何定义时间相关的离散参数?的处理方法,对大家解决问题具有一定的参考价值,需要的朋友们下面随着小编来一起学习吧!

问题描述

最近,我用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=0X>8q_in=2.5X<=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=0when 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.

这篇关于如何定义时间相关的离散参数?的文章就介绍到这了,希望我们推荐的答案对大家有所帮助,也希望大家多多支持!

08-11 15:46