2

Is there any way to model a variable which is an average value of the previous time step (t-1)'s value and the current time step (t)'s value (e.g., eq1) x(t) = (some equation), eq2) a(t) = (x(t-1) + x(t))/2) in IMODE=6 in Gekko?

The simplified code below shows what I asked above (The code itself doesn't (physically) make sense, but the code is simplified from the original to ask the question in a more succinct manner).

import numpy as np
from gekko import GEKKO

# Create GEKKO model
m = GEKKO(remote=False)
print(m.path)

m.time = np.linspace(0,5,6)

# Variables
T_air_set_prev = m.Var(value=10, name="T_air_set_prev") # previous -> previous time step (t-1)
T_m0 = m.Var(value=15, name="T_m0")
T_m_ac_t = m.Var(value=25, name="T_m_ac_t") # t -> current time step (t)
T_m_ac = m.Var(value=25, name="T_m_ac") 
Q_heating =  m.Var(value=0, name="Q_heating")

actuator = m.MV(value=6, lb=3, ub=9, name="actuator")
actuator.STATUS = 1
actuator.DMAX = 1.5

# Equations
m.Equations([T_m0 == 0.5*T_air_set_prev + 1,\
             T_m_ac_t == 0.7* T_m0 + 10,\
             T_m_ac == (T_air_set_prev + T_m_ac_t)/2,\
             Q_heating == 5*T_m_ac+20])

# This line is meant for saving the previous time step's value, but I don't know the correct way to implement this.
T_air_set_prev = m.Intermediate(T_m_ac_t)

m.Obj(0.5*Q_heating + actuator)

# Solve
m.options.IMODE = 6
m.options.DIAGLEVEL = 4
m.options.MAX_ITER = 1000
m.options.SOLVER = 3
m.solve(disp=True)
Yun Joon Jung
  • 111
  • 4

2 Answers2

2

Try using the delay function in Gekko as shown in the documentation. The number 1 is how many time steps (not time) so this should be what you need.

m.delay(T_m_ac_t,T_air_set_prev,1)

Here is the complete script.

import numpy as np
from gekko import GEKKO

# Create GEKKO model
m = GEKKO(remote=False)
#print(m.path)

m.time = np.linspace(0,5,6)

# Variables
# previous -> previous time step (t-1)
T_air_set_prev = m.Var(value=10, name="T_air_set_prev") 
T_m0 = m.Var(value=15, name="T_m0")
# t -> current time step (t)
T_m_ac_t = m.Var(value=25, name="T_m_ac_t") 
T_m_ac = m.Var(value=25, name="T_m_ac") 
Q_heating =  m.Var(value=0, name="Q_heating")

actuator = m.MV(value=6, lb=3, ub=9, name="actuator")
actuator.STATUS = 1
actuator.DMAX = 1.5

# Equations
m.Equations([T_m0 == 0.5*T_air_set_prev + 1,\
             T_m_ac_t == 0.7* T_m0 + 10,\
             T_m_ac == (T_air_set_prev + T_m_ac_t)/2,\
             Q_heating == 5*T_m_ac+20])

# Time delay
m.delay(T_m_ac_t,T_air_set_prev,1)

m.Obj(0.5*Q_heating + actuator)

# Solve
m.options.IMODE = 6
m.options.DIAGLEVEL = 0
m.options.MAX_ITER = 1000
m.options.SOLVER = 3
m.solve(disp=True)
TexasEngineer
  • 684
  • 6
  • 15
2

You might want to use the time delay function in GEKKO, m.delay. Basically, add a new gekko varaiable with m.delay() shifting one timestep, then you can use it in the equation as you usually do.

T_m_ac_t = m.Var(value=25, name="T_m_ac_t") # t -> current time step (t)
T_air_set_prev = m.Var(value=10, name="T_air_set_prev") # previous -> previous time step (t-1)

m.delay(T_m_ac_t, T_air_set_prev, 1)

Check out the example in the apmonitor page: Tiem Delay (Dead-Time)

Junho Park
  • 997
  • 4
  • 12