How to apply dead-time from FOPDT model to a Manipulated Variable in MPC gekko

Question

I'm trying to build a MPC model with previously estimated values of k, tau and theta in a FOPDT equation. I implemented this sugestion to estimate the dead-time using cspline: How to estimate theta value in FOPDT equation using gekko?.

However, I can't do the same thing with MPC because I can't use a MV as the second argument of cspline in order to apply the dead-time. If I use a gekko variable for this second argument instead, the MV's don't change because they are not present in the equation. I have two MV's, one CV and two disturbance variables. How can I apply the dead-time in this case? Thank you

import numpy as np
import time
import plotly.express as px 
from gekko import GEKKO
import json
from datetime import datetime

start_time_count = datetime.now().strftime('%Y-%m-%d %H:%M:%S')
print("start_time: ", start_time_count)

dist1 = d19jcSliceMPC['Cond_PID_SP'].values
mv1 = d19jcSliceMPC['Front_PID_PV'].values
mv2 = d19jcSliceMPC['Rear_PID_PV'].values
dist2 = d19jcSliceMPC['Pull_rolling'].values
cv1 = d19jcSliceMPC['Cond_Center_Btm_TC'].values

run_time = 3.0 * 60.0 
n = int(0.2*run_time)

Tsp1 = 1163.0 #setpoint

m = GEKKO(name='MPCbtmTC',remote=False)
m.time = np.linspace(0,n-1,n)
time_uc = m.Param(m.time)

# MV
Front = m.MV(value=mv1)
Front_ss = m.Param(value=mv1[0])
KpFront = m.Param(value=1.685312)
tauFront = m.Param(value=5.770839)
thetaFront = m.Param(value=0.114705)
t1 = np.linspace(-1,n-1,n)
ucFront = m.Var(); tm1 = m.Var(); m.Equation(tm1==time_uc-thetaFront)
m.cspline(tm1,ucFront,t1,np.array(Front),bound_x=False)

Rear = m.MV(value=mv2)
Rear_ss = m.Param(value=mv2[0])
KpRear = m.Param(value=0.1)
tauRear = m.Param(value=36.0)
thetaRear = m.Param(value=3.779397)
t2 = np.linspace(-4,n-1,n)
ucRear = m.Var(); tm2 = m.Var(); m.Equation(tm2==time_uc-thetaRear)
m.cspline(tm2,ucRear,t2,np.array(Rear),bound_x=False)

Front.STATUS = 1  # use to control temperature
Front.FSTATUS = 0 # no feedback measurement
Front.LOWER = 1160.0
Front.UPPER = 1200.0
Front.DMAX = 2.0
Front.COST = 0.0
Front.DCOST = 1.0e-4

Rear.STATUS = 1  # use to control temperature
Rear.FSTATUS = 0 # no feedback measurement
Rear.LOWER = 1180.0
Rear.UPPER = 1220.0
Rear.DMAX = 2.0
Rear.COST = 0.0
Rear.DCOST = 1.0e-4

# Parameters (disturbance)
CondSP = m.Param(value=dist1)
CondSP_ss = m.Param(value=dist1[0])
KpCondSP = m.Param(value=4.990293)
tauCondSP = m.Param(value=29.272660)
thetaCondSP = m.Param(value=2.554230)
t3 = np.linspace(-3,n-1,n)
ucCondSP = m.Var(); tm3 = m.Var(); m.Equation(tm3==time_uc-thetaCondSP)
m.cspline(tm3,ucCondSP,t3,dist1,bound_x=False)

Pull = m.Param(value=dist2)
Pull_ss = m.Param(value=dist2[0])
KpPull = m.Param(value=0.151304)
tauPull = m.Param(value=4.128567)
thetaPull = m.Param(value=0.0)
t4 = np.linspace(-0,n-1,n)
ucPull = m.Var(); tm4 = m.Var(); m.Equation(tm4==time_uc-thetaPull)
m.cspline(tm4,ucPull,t4,dist2,bound_x=False)

# Controlled variable
TC1_ss = m.Param(value=cv1[0])
TC1 = m.CV(value=TC1_ss.value)
TC1.STATUS = 1     # minimize error with setpoint range
TC1.FSTATUS = 1    # receive measurement
TC1.TR_INIT = 2    # reference trajectory
TC1.TAU = 2       # time constant for response

# Equation
m.Equation(TC1.dt()==(KpFront*(ucFront-Front_ss)-(TC1-TC1_ss))/tauFront + (KpRear*(ucRear-Rear_ss)-(TC1-TC1_ss))/tauRear+ 
                     (KpCondSP*(ucCondSP-CondSP_ss)-(TC1-TC1_ss))/tauCondSP + (KpPull*(ucPull-Pull_ss)-(TC1-TC1_ss))/tauPull)

# Global Options
m.options.IMODE   = 6 # MPC
m.options.CV_TYPE = 1 # Objective type
m.options.NODES   = 3 # Collocation nodes
m.options.SOLVER  = 1 # 1=APOPT, 3=IPOPT

TC1.MEAS = cv1[0]
# input setpoint with deadband +/- DT
DT = 0.1
TC1.SPHI = Tsp1 + DT
TC1.SPLO = Tsp1 - DT
# solve MPC
m.solve(disp=False)    
# get additional solution information
with open(m.path+'//results.json') as f:
    results = json.load(f)

finish_time_count = datetime.now().strftime('%Y-%m-%d %H:%M:%S')
print("end_time: ", finish_time_count)

df_plot = pd.DataFrame({'DateTime' : d19jcSliceMPC.time,
                        'CV (TC1)' : results['v9.bcv'],
                        'SPHI' : results['v9.sp_hi'],
                        'SPLO' : results['v9.sp_lo'],
                        'Front' : Front,
                        'Rear' : Rear})
figGekko = px.line(df_plot, 
                   x='DateTime',
                   y=['CV (TC1)','SPHI','SPLO','Front','Rear'],
                   labels={"value": "Degrees Celsius"},
                   title = "MPC")
figGekko.update_layout(legend_title_text='')
figGekko.show()

Edit:

As suggested, i changed to

ucFront = m.Var(); m.Equation(ucFront==Front) 
tm1 = m.Var(); m.Equation(tm1==time_uc-thetaFront)
m.cspline(tm1,ucFront,t1,np.array(Front),bound_x=False)

but I get this error:

Error: Exception: Access Violation
At line 359 of file ./f90/cqp.f90
Traceback: not available, compile with -ftrace=frame or -ftrace=full

Error: 'results.json' not found. Check above for additional error details

If I leave just the MV as the unshifted input I get the same error as before

TypeError: y_data must be a python list or numpy array

Answer 1

As Junho Park correctly observed, you can create another variable x such as:

MV = m.MV()  # use MV in MPC application

x  = m.Var() # use x with cspline function
m.Equation(x==MV) # connect x and MV

The cspline object complains about the MV because it wants a calculated variable instead of a type that is by default fixed and consumes a degree of freedom.