MPC con modelo ARX usando Gekko

Sep 04 2020

Estoy modelando un MPC para controlar un refrigerador y mantener la temperatura dentro de un intervalo dado mientras minimizo el costo. Estoy usando GEKKO para modelar mi algoritmo.

Escribí el siguiente código. Primero, identifiqué mi modelo usando datos de sensores de mi sistema (usé la función sysif de GEKKO). Luego construí un modelo ARX (usando la función arx en GEKKO) que se convierte en los resultados de sysid () como entrada.

Estoy tratando de escribir un algoritmo "ficticio" para probarlo localmente antes de implementarlo en un Pi.

Obtuve el siguiente error :

KeyError                                  Traceback (most recent call last)
<ipython-input-13-108148376700> in <module>
    107 #Solve the optimization problem.
    108 
--> 109 m.solve()

~/opt/anaconda3/lib/python3.8/site-packages/gekko/gekko.py in solve(self, disp, debug, GUI, **kwargs)
   2214         if timing == True:
   2215             t = time.time()
-> 2216         self.load_JSON()
   2217         if timing == True:
   2218             print('load JSON', time.time() - t)

~/opt/anaconda3/lib/python3.8/site-packages/gekko/gk_post_solve.py in load_JSON(self)
     48                             vp.__dict__[o] = dpred
     49                 else: #everything besides value, dpred and pred
---> 50                     vp.__dict__[o] = data[vp.name][o]
     51     for vp in self._variables:
     52         if vp.type != None: #(FV/MV/SV/CV) not Param or Var

KeyError: 'int_p6'

Y este es mi codigo

from gekko import GEKKO
import numpy as np
import matplotlib.pyplot as plt
m = GEKKO(remote = True)


#initialize variables

#Room Temprature:
T_external = [23,23,23,23,23.5,23.5,23.4,23.5,23.9,23.7,\
              23,23.9,23.9,23.4,23.9,24,23.6,23.7,23.8,\
              23,23,23,23,23]

# Temprature Lower Limit:
temp_low = 10*np.ones(24)

# Temprature Upper Limit:
temp_upper = 12*np.ones(24)

#Hourly Energy prices:
TOU_v = [39.09,34.93,38.39,40.46,40.57,43.93,25,11,9,24,51.28,45.22,45.72,\
            36,35.03,10,12,13,32.81,42.55,8,29.58,29.52,29.52]

###########################################
#System Identification:

#Time 
t = np.linspace(0,10,117)
#State of the Fridge
ud = np.append(np.zeros(78) ,np.ones(39),0)
#Temprature Data
y = [14.600000000000001,14.600000000000001,14.700000000000001,14.700000000000001,14.700000000000001,\
     14.700000000000001,14.700000000000001,14.700000000000001,14.700000000000001,14.700000000000001,\
     14.700000000000001,14.700000000000001,14.700000000000001,14.8,14.8,14.8,14.8,14.8,14.8,14.8,14.8,\
    14.8,14.8,14.9,14.9,14.9,14.9,14.9,14.9,14.9,15,15,15,15,15,15,15,15,15,15,15,15,15.100000000000001,\
    15.100000000000001,15.100000000000001,15.100000000000001,15.100000000000001,15.100000000000001,\
    15.100000000000001,15.100000000000001,15.100000000000001,15.100000000000001,15.100000000000001,\
    15.100000000000001,15.100000000000001,15.100000000000001,15.100000000000001,15.100000000000001,\
    15.100000000000001,15.100000000000001,15.100000000000001,15.100000000000001,15.100000000000001,\
    15.100000000000001,15.100000000000001,15.100000000000001,15.100000000000001,15.100000000000001,\
    15.100000000000001,15.100000000000001,15.100000000000001,15.100000000000001,15.100000000000001,\
    15.100000000000001,15.100000000000001,15.100000000000001,15.100000000000001,15.100000000000001,
    15,15,15,15,15,15,15,15,15,15,14.9,14.9,14.9,14.9,14.8,14.9,14.8,14.8,14.8,14.8,14.8,14.8,\
    14.8,14.700000000000001,14.8,14.700000000000001,14.700000000000001,14.700000000000001,\
    14.700000000000001,14.700000000000001,14.700000000000001,14.700000000000001,\
    14.700000000000001,14.600000000000001,14.600000000000001,14.600000000000001,\
    14.600000000000001,14.600000000000001,14.60]

na = 1 # output coefficients
nb = 1 # input coefficients
print('Identification')
yp,p,K = m.sysid(t,ud,y,na,nb,objf=10000,scale=False,diaglevel=1)
#create control ARX model:

y = m.Array(m.CV,1)
uc = m.Array(m.MV,1)
m.arx(p,y,uc)
# rename CVs
T= y[0]

# rename MVs
uc = uc[0]

# steady state initialization
m.options.IMODE = 1
m.solve(disp=True)

###########################################

#Parameter
P = m.Param(value =100) #power
TL = m.Param(value=temp_low) 
TH = m.Param(value=temp_upper)
c = m.Param(value=TOU_v)
# Manipilated variable:

u = m.MV(lb=0, ub=1, integer=True)
u.STATUS = 1  # allow optimizer to change the variable to attein the optimum.

# Controlled Variable (Affected with changes in the manipulated variable)

T = m.CV(value=11) # Temprature will start at 11.

# Soft constraints on temprature.

eH = m.CV(value=0)
eL = m.CV(value=0)

eH.SPHI=0       #Set point high for linear error model.
eH.WSPHI=100    #Objective function weight on upper set point for linear error model.
eH.WSPLO=0      # Objective function weight on lower set point for linear error model
eH.STATUS =1    # eH : Error is considered in the objective function.
eL.SPLO=0
eL.WSPHI=0
eL.WSPLO=100 
eL.STATUS = 1   
#Linear error (Deviation from the limits)
m.Equations([eH==T-TH,eL==T-TL])

#Objective : minimize the costs.

m.Minimize(c*P*u)

#Optimizer Options.

m.options.IMODE = 6 # MPC mode in Gekko.
m.options.NODES = 2  # Collocation nodes.
m.options.SOLVER = 1 # APOT solver for mixed integer linear programming.
m.time = np.linspace(0,23,24)

#Solve the optimization problem.

m.solve()

Respuestas

2 JohnHedengren Sep 04 2020 at 03:41

El problema es con:

T = m.CV(value=11) # Temperature will start at 11.

Está redefiniendo la Tvariable, pero almacena ambas internamente. Si necesita reinicializar, 11use T.value=11. Además, agregué las variables eHy eLantes de la inicialización del estado estable. Aquí hay un script completo que se ejecuta correctamente.