Gekko to simulate process with disturbances

Question

I'm learning GEKKO and I was wondering how to simulate a process and then control with some disturbances variables involved.

Let's pick a simple heat exchanger for example, where I'm trying to heat up a water stream from 100ºC to 110ºC using a heat source from a steam of 6,5 bar.

Let's assume the only heat source comes from the condensation of steam, so our mass balance would be:

Q = m_s * lamb_s = m_w * Cp * d(T_out - T-in)dt

Where m_s: steam mass flow lamb_s: steam latent heat m_w: water mass flow Cp: water specific heat T_out: Water temperature out T_in: water temperature in

Now, I want to control T_out with m_s with a constant setpoint for T_out of 110ºC. But I want to be able to step m_w and T_in. I get that it's easy to step m_w because it is indepent of time, but how do I step T_in if in the equation it is time relative (dT_in/dt) ? If I'm stepping T_in, I cant cut it from the equation, right?

THis is what I've tried this, where I thought that dT_in/dt would only not be 0 when I'm stepping it.:

Tia[300] = -10 +273.15 Ti_in = m.Param(value=Tia) # K

m.Equation(T_out.dt() == (1.0/(m_wCp))(m_s * lamb_s)) + Ti_in)

BUt the response was a step in T_out.

I also tried writing in the equation Ti_in as T_in.dt(), but I got AttributeError: dt, I guess it's because I'm defining it as a parameter?

Answer 1

The equation for an energy balance is:

m.Equation(mass*Cp*T_out.dt() == m_s*lambda_s 
                                 + mdot*Cp*(T_in - T_out))

where mass is the fluid mass in the control volume and T_in, m_s are independent inputs to the model. T_in can be defined as an input disturbance and m_s as the manipulated variable. Here is an example script that demonstrates this in a Model Predictive Control application with a step in the T_in disturbance.

Here is the complete script:

from gekko import GEKKO
import numpy as np

# Initialize Model
m = GEKKO(remote=False)

# Time points
n = 101
m.time = np.linspace(0, 10, n)

# Parameters
lambda_s = m.Param(value=2257)    # Steam latent heat (kJ/kg)
Cp = m.Param(value=4.18)          # Water specific heat (kJ/kg.K)
mdot = m.Param(value=1)           # Water mass flow (kg/s)
mass = m.Param(value=0.5)         # Mass of vessel fluid (kg)

# Disturbance
d = np.ones(n)*(100+273.15); d[50:]=(95+273.15)
T_in = m.Param(d)  # Water temperature in (K)

# Controlled Variable
# Water temperature out (K)
T_out = m.CV(value=100+273.15,name='tout')
T_out.SPHI = 110 + 273.15 + 0.5  # Setpoint High
T_out.SPLO = 110 + 273.15 - 0.5   # Setpoint Low
T_out.TAU = 1
T_out.TR_INIT = 1
T_out.STATUS = 1                  # Enable control

# Manipulated Variable
m_s = m.MV(value=0, lb=0, ub=5)   # Steam mass flow (kg/s)
m_s.DCOST = 1e-3
m_s.STATUS = 1                    # Allow optimizer to change m_s
m_s.FSTATUS = 0                   # No MV measurement feedback

# Equation
m.Equation(mass*Cp*T_out.dt() == m_s*lambda_s 
                                 + mdot*Cp*(T_in - T_out))

# Solve MPC
m.options.IMODE = 6   # Control mode
m.options.CV_TYPE = 1 # l1-norm objective
m.options.SOLVER = 3  # Solver
m.solve(disp=True)

# Plot results
# get additional solution information
import json
with open(m.path+'//results.json') as f:
    results = json.load(f)

import matplotlib.pyplot as plt
plt.figure()
plt.subplot(3,1,1)
plt.plot(m.time,T_in.value,color='orange',
         linestyle='-',label='DV Step')
plt.ylabel(r'$T_{in}$'); plt.legend(); plt.grid()
plt.subplot(3,1,2)
plt.plot(m.time,m_s.value,'b-',label='MV Optimized')
plt.ylabel(r'$m_s$'); plt.legend(); plt.grid()
plt.subplot(3,1,3)
plt.plot(m.time,results['tout.tr_hi'],'k-.',label='SPHI Traj')
plt.plot(m.time,T_out.value,'r--',label='CV Response')
plt.plot(m.time,results['tout.tr_lo'],'k:',label='SPLO Traj')
plt.ylabel(r'$T_{out}$'); plt.xlabel('Time'); plt.grid()
plt.legend(loc='best'); plt.tight_layout()
plt.savefig('mpc.png',dpi=300); plt.show()