▷ Sintonía del Controlador MPC (MPC Tuning)

La sintonía de controladores predictivos basados en modelos (Model Predictive Control, MPC) es una tarea crucial en ingeniería de control avanzada, ya que determina el rendimiento dinámico, la estabilidad y la robustez del sistema controlado. Este artículo profundiza en los principios, los parámetros de ajuste y los métodos avanzados utilizados para la optimización de MPC, empleando la formulación matemática estándar y ejemplos prácticos relevantes.

Antes de continuar te invito a que te inscribas al canal

Contenidos mostrar

Tuning of Model Predictive Controllers Based on Hybrid Optimization

Antes de comenzar a entender sobre la sintonia del controlador Predictivo, te dejo el trabajo empleado como referencia para la elaboración de este post el cual adicionalmente propone un método de sintonia por medio de optimización hibrida para controladores predictivos basados en modelo MPC y NMPC.

Descargar Articulo

Los códigos de implementación junto con diversos ejemplos puedes descargarlos del siguiente repositorio los cuales están adaptados para emplear el propio MPC Toolbox de Matlab, basta con entrar a la carpeta MPC-Tuning.

Descargar Códigos

Principios del Control Predictivo Basado en Modelos

El control predictivo basado en modelos (MPC) es una estrategia avanzada ampliamente utilizada en procesos industriales y sistemas complejos. Su principal fortaleza radica en su capacidad para manejar restricciones explícitas, optimizar el desempeño dinámico y operar con modelos multivariables. A continuación, se describen los principios fundamentales de MPC, incluyendo su formulación matemática, naturaleza predictiva y enfoque de optimización.

Formulación del Problema de Control Predictivo

El MPC utiliza un modelo dinámico explícito del sistema para predecir su comportamiento futuro. En cada intervalo de muestreo, el controlador resuelve un problema de optimización basado en una función de costo, buscando minimizar el error entre las salidas predichas y las referencias deseadas, así como las variaciones en las entradas manipuladas.

La formulación clásica del problema de optimización en MPC es la siguiente:

\underset{\mathbf{ u, \epsilon_k}}{\rm min}\ \Bigg\{ {J}=\sum_{i=1}^{n_y}\sum_{n=1}^{p}\left \{ {Q}_{ii} \left[ \hat{{y}}_i(k+n|k)-{r}_i(k+n) \right]^2 \right \}+

\sum_{j=1}^{n_u}\sum_{n=1}^{m_j}\left \{{W}_{jj}\left[ {\Delta u}_j(k+n-1)\right]^2 \right \}+ \rho \epsilon_k^2\Bigg\},

{\rm subject\ to}\\\\ x(k+n)=f(x(k+n-1),\ u(k+n),\ \boldsymbol{\mu} ),\ n=1,\ldots,p\\ g(y(k+n),\ x(k+n),\ u(k+n),\ \boldsymbol{\mu})=\mathbf{0},\ n=1,\ldots,p\\ x(k) = \hat{x}(k)\\ {u_{j_{min}}}\leq {u}_j(k+n-1)\leq {u_{j_{max}}},\ j=1,\ldots,n_u\ {\rm and}\ n = 1, \ldots, m_j\\ {\Delta u_{j_{min}}}\leq {u}_j(k+n-1)-{u}_j(k+n-2)\leq \Delta {u_{j_{max}}},\\ j=1,\ldots,n_u\ {\rm and}\ n = 1, \ldots, m_j\\ {y_{i_{min}}}\leq \hat{{y}}_i(k+n)\leq {y_{i_{max}}},\ i=1,\ldots,n_y\ {\rm and}\ n = 1, \ldots, p,

Donde:

J: Función de costo que se optimiza en cada intervalo de muestreo.
$\mathbf{u}$ : Vector de entradas manipuladas.
$\mathbf{y}$ : Vector de salidas controladas.
$\mathbf{x}$ : Vector de variables de estado.
$\hat{y}_i(k+n|k)$ : Valor predicho de la salida i-ésima en el paso n.
$r_i(k+n)$ : Trayectoria de referencia deseada para la salida i-ésima en el paso n.
$\mathbf{Q}$ : Matriz diagonal de pesos para los errores de las salidas.
$\mathbf{W}$ : Matriz diagonal de pesos para las variaciones en las entradas manipuladas.
$\epsilon_k$ : Variable de relajación para restricciones suaves.
$\rho$ : Penalización asociada a $\epsilon_k$ .

Principios Predictivos

El núcleo del MPC radica en su capacidad de predecir el comportamiento futuro del sistema utilizando un modelo dinámico, típicamente representado por ecuaciones de estado o de transferencia. A partir de las condiciones iniciales y las entradas actuales, el modelo genera una predicción de las salidas del sistema para un horizonte futuro (p).

Modelo de Estado

El modelo dinámico de estado se representa como:

x(k+1)=Ax(k)+Bu(k)\\y(k)=Cx(k)+Du(k)

donde:

A,B,C,D: Matrices del modelo dinámico del sistema.
x(k): Estado del sistema en el instante k.
u(k): Entrada manipulada en el instante k.
y(k): Salida controlada en el instante k.

Este modelo puede extenderse para incluir incertidumbres y no linealidades mediante formulaciones más avanzadas, como MPC no lineal (NMPC) o robusto.

Parámetros Clave para la Sintonía del MPC

La sintonía del MPC implica ajustar parámetros fundamentales que determinan el desempeño dinámico del sistema, la estabilidad, la robustez y el cumplimiento de restricciones. Este punto aborda con detalle los cuatro parámetros esenciales: el horizonte de predicción (p), el horizonte de control ( $m_j$ ), las matrices de pesos ( $Q_{ii}$ y $W_{jj}$ ), y el periodo de muestreo ( $T_s$ ).

Horizonte de Predicción (p)

El horizonte de predicción (p) define cuántos pasos hacia adelante el controlador considera para predecir el comportamiento del sistema. Este parámetro es fundamental para capturar la dinámica del proceso y optimizar el rendimiento del sistema.

Impacto del Horizonte de Predicción

Horizonte largo:
Proporciona una mayor capacidad para anticiparse a los efectos dinámicos y a las restricciones futuras, especialmente útil para sistemas lentos o con dinámicas complejas. Sin embargo, incrementa el costo computacional, ya que el problema de optimización se resuelve sobre un espacio de búsqueda más amplio.
Horizonte corto:
Es más adecuado para sistemas rápidos o con modelos altamente inciertos, donde las predicciones a largo plazo pueden ser menos confiables. Además, reduce significativamente la carga computacional.

Horizonte de Control ( $m_j$ )

El horizonte de control ( $m_j$ ) representa el número de movimientos futuros de la entrada manipulada que se optimizan en cada iteración. Este parámetro influye directamente en la flexibilidad del controlador para ajustar las señales de control.

Impacto del Horizonte de Control

$m_j = 1$ :
Un horizonte mínimo implica que el controlador realiza ajustes inmediatos y directos, lo que reduce el costo computacional, pero puede limitar la capacidad del sistema para manejar restricciones y dinámicas complejas.
$m_j > 1$ :
Incrementa los grados de libertad del controlador, permitiendo optimizaciones más sofisticadas, pero a costa de un mayor esfuerzo computacional y, en algunos casos, mayor sensibilidad a perturbaciones.

Generalmente, se establece que: $m_j \leq p$ ,

para garantizar que las decisiones de control se alineen con las predicciones. Sin embargo, en sistemas altamente dinámicos, un $m_j$ mayor puede ser beneficioso para lograr respuestas más suaves.

Matrices de Pesos ( $Q_{ii}$ y $W_{jj}$ )

Las matrices de pesos son fundamentales para priorizar ciertos objetivos del controlador y penalizar comportamientos no deseados.

La matriz diagonal $\mathbf{Q}$ define la penalización por error entre las salidas predichas y las referencias. Los pesos se asignan según la importancia relativa de cada variable controlada:

Pesos altos:
Priorizan variables críticas, como la presión o temperatura en procesos químicos sensibles.
Pesos bajos:
Reducen la prioridad de ciertas salidas, permitiendo al sistema enfocarse en otros objetivos.

La matriz diagonal $\mathbf{W}$ penaliza los cambios en las entradas manipuladas para evitar movimientos bruscos:

Pesos altos:
Suavizan las variaciones en las señales de control, reduciendo el desgaste de actuadores y aumentando la estabilidad.
Pesos bajos:
Permiten cambios más agresivos, útiles en sistemas que requieren respuestas rápidas.

Periodo de Muestreo ( $T_s$ )

El periodo de muestreo ( $T_s$ ) determina la frecuencia con la que se actualizan las predicciones y decisiones del controlador. Un $T_s$ adecuado es crucial para evitar fenómenos de aliasing y reducir el costo computacional.

El periodo de muestreo debe ser:

Lo suficientemente corto:
Para capturar las dinámicas rápidas del sistema.
No demasiado corto:
Para evitar aumentar innecesariamente la carga computacional y el ruido en las señales de medición.

Estrategia MPCTuning: Optimización Híbrida para la Sintonía del MPC

La estrategia MPCTuning es un método híbrido de optimización diseñado para la sintonía de controladores predictivos basados en modelos (MPC) en sistemas multi-input multi-output (MIMO). Su enfoque combina el método de logro de objetivos (Goal Attainment Method, GAM) y la búsqueda de vecindario variable (Variable Neighborhood Search, VNS), lo que permite encontrar los parámetros óptimos del controlador minimizando el error de seguimiento y garantizando un bajo costo computacional.

Metodología de Optimización

La sintonía del MPC se lleva a cabo en dos etapas:

Optimización de las matrices de pesos del MPC con el método GAM:

Se ajustan las matrices de pesos Q y W en la función objetivo del controlador.
Se minimiza el error cuadrático entre la respuesta del modelo interno del MPC en lazo cerrado y una trayectoria de referencia predefinida.
Se plantea un problema de optimización multiobjetivo, buscando una solución de compromiso en la frontera de Pareto.

Optimización de los horizontes de predicción y control mediante VNS:

Se determinan los valores óptimos de los horizontes de predicción (p) y control (m).
Se minimiza el error cuadrático entre la respuesta en lazo cerrado y una trayectoria óptima.
Se busca un equilibrio entre bajo costo computacional y alto desempeño.

Método de Logro de Objetivos (GAM)

El método GAM permite resolver problemas de optimización con múltiples objetivos en conflicto. En la sintonía del MPC, se define un conjunto de objetivos ideales:

\mathbf{f^*} = [f^*_1, f^*_2, … , f^*_n]

Estos corresponden a una solución utópica inalcanzable. El problema se reformula introduciendo pesos $\omega_i$ y variables de holgura para encontrar una solución realista:

\underset{\mathbf{x},\gamma}{\rm min}\ \gamma,\\ {\rm subject\ to}\\ \mathbf{h}(\mathbf{x})=0\\ \mathbf{g}(\mathbf{x})\leq 0\\ f_i(\mathbf{x})-\omega_i\gamma-f_i^*\leq 0,\ i=1,\cdots,n\\ \mathbf{L_B}\leq \mathbf{x} \leq \mathbf{U_B}

[
\min_{\mathbf{x}, \gamma} \quad \gamma, \quad \text{sujeto a} \quad f_i(x) – \omega_i \gamma – f^*_i \leq 0, \quad i = 1, … , n
]

Este enfoque permite obtener un balance entre desempeño y robustez, ajustando dinámicamente las matrices Q y W para priorizar ciertos objetivos del controlador.

Búsqueda de Vecindario Variable (VNS)

El VNS es un algoritmo metaheurístico que explora diferentes soluciones dentro de un vecindario y evita caer en óptimos locales. En la sintonía del MPC, el algoritmo:

Explora diferentes configuraciones de horizontes de predicción y control.
Evalúa la calidad de cada configuración en función del error de seguimiento y la complejidad computacional.
Realiza cambios progresivos en los valores de p y m para encontrar la mejor combinación.

La formulación del problema de optimización con VNS es:

\underset{p,\mathbf{m}}{\rm min}\ \bigg\{ f_v(\mathbf{x}_{dv})= \sum_{i=1}^{n_y}\left[ \sum_{k=1}^{\phi}\left( \left \{ {y}_i(k)-{y}_{o_i}(k|1) \right \}^2 + \left \{ {y}_i^{R}(k)-{y}_i(k) \right \}^2\right)\right]\\ +p+ \sum_{j=1}^{n_u}\sum_{k=1}^{m_j} \left[ \dfrac{|u_{o_j}(1|1)|}{|u_{o_j}(k+1|1)-u_{o_j}(k|1)|} \right]^2 \bigg \}\\ {\rm subject\ to}\\ Eq. MPC\\ u_{o_j}(k+1|1)\neq u_{o_j}(k|1)\\ \mathbf{m}< p,

[
\min_{p, m} \sum_{i=1}^{n_y} \sum_{k=1}^{\phi} \left( (y_i(k) – y_o_i(k|1))^2 + (y_R_i(k) – y_i(k))^2 \right) + p + \sum_{j=1}^{n_u} \sum_{k=1}^{m_j} |u_o_j(1|1) – u_o_j(k|1)|
]

Donde:

(y(k)) es la salida en lazo cerrado del controlador.
(y_o(k|1)) es la primera predicción del MPC en el instante inicial.
(y_R(k)) es la trayectoria de referencia predefinida.
(p) y (m) son los horizontes de predicción y control, respectivamente.

Este método permite encontrar un equilibrio entre precisión y costo computacional, evitando configuraciones excesivamente complejas que no aporten mejoras significativas al desempeño del sistema.

Evaluación de Desempeño y Robustez

Para garantizar que el método MPCTuning proporcione controladores robustos, se incorporan estrategias como:

Uso de formulaciones robustas en el MPC y en el MPCT, considerando incertidumbres en el modelo.
Estimación de parámetros en escenarios de peor caso, utilizando una familia de modelos representativos.
Minimización de errores en la trayectoria de referencia, asegurando un desempeño consistente ante variaciones en el proceso.

El MPCTuning ha sido probado en procesos industriales de referencia, demostrando mejoras en el desempeño y reducción del esfuerzo computacional en comparación con métodos de sintonía tradicionales.

A continuación, se van a explorar varios ejemplos empleando la metodología MPC Tuning, ilustrando la aplicación práctica de diferentes configuraciones de control predictivo multivariable.

Subsistema MIMO 3×3 de la Shell Heavy Oil Fractionator

Un subsistema MIMO $3 \times 3$ del sistema Shell Heavy Oil Fractionator (HOF) propuesto en Maciejowski se sintonizó para mostrar el desempeño del algoritmo $\text{MPCT}$ . Para este caso, se eligió una formulación lineal de MPC, conocida como Generalized Predictive Control (GPC), aunque cualquier otra formulación lineal también sería aplicable.

La misma consigna (o setpoint) se ha configurado para el algoritmo GAM y para la planta, a fin de contemplar la interacción interna entre las variables del sistema MIMO. Las restricciones de entrada y de incrementos de entrada para el controlador MPC se definen en LaTeX como sigue:

\mathbf{u_{min}} = [-0.5, -0.5, -0.5]\\ \mathbf{u_{max}} = [0.5, 0.5, 0.5]\\ \mathbf{\Delta u_{min}} = [-0.05, -0.05, -0.05]\\ \mathbf{\Delta u_{max}} = [0.05, 0.05, 0.05]

Los setpoints se cambian en los instantes de tiempo:

70 min: $\mathbf{w} = [0.2, 0.2, 0.2]$
315 min: $\mathbf{w} = [0.0, 0.4, 0.1]$
800 min: $\mathbf{w} = [0.1, 0.3, 0.0]$
1600 min: $\mathbf{w} = [0.0, 0.0, 0.0]$

Además, se introducen perturbaciones desconocidas tipo escalón de intensidad -0.05 en la entrada $u_1$ (desde 1100 min hasta 1120 min), y de intensidad 0.1 en la entrada $u_2$ (desde 1400 min hasta 1420 min).

Las tres entradas del sistema, $u_1, u_2, u_3$ , son:

Caudal de extracción superior
Caudal de extracción lateral
Carga térmica en el rehervidor de la parte inferior

Mientras que las salidas controladas, $y_1, y_2, y_3$ , representan:

Composición del punto final superior
Composición del punto final lateral
Temperatura en el rehervidor inferior

El modelo dinámico de este sistema se describe mediante la siguiente matriz de funciones de transferencia:

{\mathbf{G}}(s)= \left[ {\begin{matrix} {\dfrac{{ 4.05+2.11\epsilon_1}}{{50s+1}}e^{ -27 s}}&{\dfrac{{ 1.77+0.39\epsilon_2}}{{60s+1}}e^{ -28s}} & {\dfrac{{ 5.88+0.59\epsilon_3}}{{50s+1}}e^{ -27 s}}\\\\ {\dfrac{{ 5.396+3.29\epsilon_1}}{{50s + 1}}e^{ - 18s}}&{\dfrac{{ 5.72+0.57\epsilon_2}}{{60s + 1}}e^{ - 14s}}&{\dfrac{{ 6.90+0.89\epsilon_3}}{{40s+1}}e^{ -15 s}}\\\\ {\dfrac{{ 4.38+3.11\epsilon_1}}{{33s+1}}e^{ -20 s}}&{\dfrac{{ 4.42+0.73\epsilon_2}}{{44s+1}}e^{ -22 s}}&{\dfrac{{ 7.20+1.33\epsilon_3}}{{19s+1}}} \end{matrix}} \right],

donde cada $\epsilon_i$ representa la incertidumbre en la ganancia del modelo y las constantes de tiempo están dadas en minutos.

Columna de Separación Binária de Wood and Berry

El sistema estudiado es una columna de destilación para la separación de una mezcla binaria de metanol y agua. La columna de destilación opera con dos entradas de control:

$u_1$ : la tasa de reflujo en la parte superior.
$u_2$ : la tasa de flujo de vapor hacia el rehervidor en la parte inferior.

Los objetivos del control son mantener las fracciones molares de metanol en el destilado ( $y_1$ ) y en el residuo ( $y_2$ ) en valores deseados.

El modelo de Wood & Berry está representado por la siguiente matriz de funciones de transferencia:

G(s) = \begin{bmatrix} \dfrac{12.8 e^{-s}}{16.7s + 1} & \dfrac{-18.9 e^{-3s}}{21.0s + 1} \\\\ \dfrac{6.6 e^{-7s}}{10.9s + 1} & \dfrac{-19.4 e^{-3s}}{14.4s + 1} \end{bmatrix}

Donde:

$y_1$ : Fracción molar de metanol en el tope de la columna (destilado).
$y_2$ : Fracción molar de metanol en el fondo de la columna (residuo).
$u_1$ : Tasa de reflujo en la parte superior.
$u_2$ : Flujo de vapor hacia el rehervidor.

Adicionalmente, se modela la perturbación medida como:

G_d(s) = \begin{bmatrix} \dfrac{3.8 e^{-8.1s}}{14.9s + 1} \\\\ \dfrac{4.9 e^{-3.4s}}{13.2s + 1} \end{bmatrix} d

Donde $d$ representa la variación en la tasa de alimentación de la columna, la cual es un disturbio medido que puede afectar las fracciones molares del destilado y del residuo.

El objetivo del sistema de control es regular $y_1$ y $y_2$ en los valores deseados mediante el ajuste de $u_1$ y $u_2$ y rechazar perturbaciones en la tasa de alimentación de la columna ( $d$ ).

Para lograr estos objetivos, sintonice adecuadamente el controlador MPC empleando el algoritmo MPCT (MPC tuning)

Reactor de Van de Vusse (NMPC)

A continuación, se expone un tercer ejemplo que ilustra el comportamiento del controlador MPC usando la sintonía propuesta, aplicado a un reactor continuo de tanque agitado (CSTR) no lineal que lleva a cabo las reacciones de Van de Vusse.

En este caso, se selecciona una formulación no lineal de MPC (NMPC) para controlar el sistema. El reactor de Van de Vusse ( $A\overset{k_1}{\rightarrow} B \overset{k_2}{\rightarrow} C$ y $2A\overset{k_3}{\rightarrow} D$ ) comprende dos reacciones del reactivo A hacia el producto deseado B y hacia los subproductos no deseados C y D.

Los parámetros de velocidad de reacción $k_i$ , $i = 1,2,3$ , dependen de la temperatura $T$ , y se representan mediante la ecuación de Arrhenius:

k_i(T)=k_{i0}\; exp\left(\dfrac{-E_i/R}{T(^\circ C)+273.15}\right),

donde $E_i$ , $i = 1,2,3$ , son las energías de activación de las tres reacciones y $R$ es la constante universal de los gases.

El proceso se describe a través de un modelo no adiabático representado por los siguientes balances de masa y energía en el reactor:

\dfrac{dC_A}{dt}=\dfrac{F}{V}(C_{Af}-C_A)-k_1(T)C_A-k_3(T)C_A^2, \\ \dfrac{dC_B}{dt}=-\dfrac{F}{V}C_B+k_1(T)C_A-k_2(T)C_B, \\ \dfrac{dT}{dt}=\dfrac{1}{\rho C_p}\left[k_1(T)C_A(-\Delta H_{RAB})+k_2(T)C_B(-\Delta H_{RBC})+k_3(T)C_A^2(-\Delta H_{RAD})\right]+\dfrac{F}{V}(T_0-T)+\dfrac{K_w A_R}{\rho C_p V}(T_k-T),

donde la concentración de $A$ en el reactor y en la alimentación son $C_A$ y $C_{Af}$ , respectivamente; $C_B$ es la concentración del producto deseado B; las entradas manipuladas son la razón de dilución $F/V$ y la temperatura de la camisa del reactor $T_k$ ; $V$ es el volumen constante del reactor; $T_0$ es la temperatura de alimentación; $\rho$ es la densidad del líquido; $C_p$ es la capacidad calorífica; $Q = K_w A_R (T - T_k)$ es el calor transferido desde el reactor hacia la camisa, donde $K_w$ es el coeficiente de transferencia de calor y $A_R$ el área de transferencia de calor; los calores de reacción se denotan por $(-\Delta H_{RAB}), (-\Delta H_{RBC}), (-\Delta H_{RAD})$ .

Los valores de los parámetros del sistema fueron obtenidos de Trierweiler y se describen en la entrada del modelo no lineal de Van de Vusse.

Las dos entradas del sistema, $u_1$ y $u_2$ , son la razón de dilución $F/V$ y la temperatura de la camisa del reactor $T_k$ , respectivamente. Las dos salidas controladas, $y_1$ y $y_2$ , corresponden a la concentración en $C_B$ y a la temperatura del reactor $T$ . Las parejas entrada/salida del reactor de Van de Vusse se establecen como $u_1:y_1$ y $u_2:y_2$ . El periodo de muestreo del NMPC se fija en $T_s=0.05$ horas. El horizonte de predicción se establece en $p=32$ (5 bits), el horizonte de control se establece como $\mathbf{m}=[7,7]$ (3 bits cada uno) y el horizonte de sintonía se fija en $\phi=60$ . Se consideran dos casos para la trayectoria de referencia, gobernada por sistemas de primer orden, donde las ganancias estáticas son $\mathbf{k}^{R}=[1,1]$ y las dos constantes de tiempo son $\boldsymbol{\tau_a}^{R}=[0.05,0.0875]$ horas y $\boldsymbol{\tau_c}^{R}=[0.3,0.4]$ horas para las dinámicas agresiva y conservadora, respectivamente.

Eso es todo por la entrada del dia de hoy, espero les haya gustado y hayan aprendido algo nuevo. Si te ha servido el contenido de esta entrada, de los videos y los códigos de implementación y deseas apoyar mi trabajo invitandome a un café super barato, puedes hacerlo en el siguiente link:

👉 Invitar a Sergio a un Café ☕️

Que esten muy bien, nos vemos en la siguiente entrada.

Sergio A. Castaño Giraldo

Mi nombre es Sergio Andres Castaño Giraldo, y en este sitio web voy a compartir una de las cosas que mas me gusta en la vida y es sobre la Ingeniería de Control y Automatización. El sitio web estará en constante crecimiento, voy a ir publicando material sobre el asunto desde temas básicos hasta temas un poco más complejos. Suscríbete al sitio web, dale me gusta a la página en Facebook y únete al canal de youtube. Espero de corazón que la información que comparto en este sitio, te pueda ser de utilidad. Y nuevamente te doy las gracias y la bienvenida a control automático educación.