Leveraging Input Convex Neural Networks (ICNNs), ICNN-based Model Predictive Control (MPC) successfully attains globally optimal solutions by upholding convexity within the MPC framework. However, current ICNN architectures encounter the issue of vanishing gradients, which limits their ability to serve as deep neural networks for complex tasks. Additionally, the current neural network-based MPC, including conventional neural network-based MPC and ICNN-based MPC, faces slower convergence speed when compared to MPC based on first-principles models. In this study, we leverage the principles of ICNNs to propose a novel Input Convex LSTM for Lyapunov-based MPC, with the specific goal of reducing convergence time and mitigating the vanishing gradient problem while ensuring closed-loop stability. From a simulation study of a nonlinear chemical reactor, we observed a mitigation of vanishing gradient problem and a reduction in convergence time, with a percentage decrease of 46.7%, 31.3%, and 20.2% compared to baseline plain RNN, plain LSTM, and Input Convex Recurrent Neural Network, respectively [Source code is available at <https://github.com/killingbear999/ICLSTM>].

Model Predictive Control, Deep Learning, Input Convex Neural Network, Computational Efficiency, Nonlinear Processes

§ INTRODUCTION

Neural network-based Model Predictive Control (MPC) has found application in various domains, such as machine profile direction control in paper manufacturing [Lanzetti et al., 2019], batch crystallization process [Zheng et al., 2022, Zheng et al., 2022], robot manipulator control [Nubert et al., 2020], regulation of Heating, Ventilation and Air-conditioning (HVAC) systems [Afram et al., 2017, Ellis and Chinde, 2020], approximation of the hybrid neuroprosthesis system [Bao et al., 2017], thin-film decomposition of quantum dot [Sitapure and Kwon, 2022]. Conventional MPC relies on the development of first-principles models, a process known to be resource-intensive. In today's big data era, data-driven deep learning approaches have emerged as viable substitutes for these first-principles models within MPC formulations. This advancement paves the way for the practical application of MPC in various industries, significantly enhancing its commercial viability.

However, employing a conventional neural network to capture system dynamics within MPC can introduce non-convexity due to the inherent non-convex nature of neural networks, often resulting in suboptimal local solutions. To ensure the attainment of a globally optimal solution and maintain the convexity of neural network-based MPC, Input Convex Neural Networks (ICNNs) have emerged as a powerful tool in the realm of optimization and control. The concept of ICNNs was initially introduced by Amos et al., 2017, and was subsequently explored and expanded by Chen et al., 2018 and Bünning et al., 2021. ICNNs have been applied to several neural network-based optimization problems, such as the Van de Vusse reactor [Yang and Bequette, 2021], optimal transport mapping [Makkuva et al., 2020], voltage regulation [Chen et al., 2020, Chen et al., 2020], molecular discovery [Alvarez-Melis et al., 2021], and DC optimal power flow [Zhang et al., 2021].

ICNN-based MPC presents a promising solution to address the challenges encountered in conventional neural network-based MPC. However, we observe that current ICNN architectures suffer from the vanishing gradient problem. This problem limits their ability to learn complex dynamics, primarily because they struggle to be trained as deep neural networks.
Moreover, in all previous research on neural network-based MPC, including ICNN-based MPC, the performance assessment has primarily centered on comparing state trajectories with established baselines and MPC with first-principles models, and ensuring the final convergence. However, in our previous works [Wu et al., 2019, Wu et al., 2019], we noted that RNN-based MPC exhibited significantly slower convergence speed compared to the MPC based on the first-principles model. In a real-life MPC for a chemical plant, numerous operations require real-time or near-real-time control to maintain product quality, safety, and operational efficiency. Swift decision-making is pivotal for safety in chemical processes, as delays in addressing reactant changes can result in undesired reactions or unsafe conditions. Rapid decision-making extends its benefits to optimizing the utilization of raw materials, energy, and other resources in other industries as well, ultimately yielding cost savings and reducing the environmental footprint. In summary, for a neural network-based MPC, the convergence runtime is a critical parameter to safeguard product quality, safety, efficiency, and resource utilization, with profound implications for both operational and environmental sustainability.

Therefore, in this study, we build upon the groundwork laid by Amos et al., 2017, Chen et al., 2018, and Bünning et al., 2021 to introduce a novel member to the family of ICNNs, namely Input Convex LSTM (ICLSTM). This extension aims to bolster the overall performance of neural network-based MPC, particularly in terms of convergence runtime, and addressing the vanishing gradient problem encountered by the current ICNNs. Additionally, we incorporate a Lyapunov-based constraint into the MPC formulation, addressing a theoretical closed-loop stability aspect that previous studies have often overlooked.

Section <ref> introduces the notation and the nonlinear Lyapunov-based Model Predictive Control (LMPC) design. Section <ref> provides a comprehensive overview of various ICNN variants and proposes a novel ICLSTM architecture, along with the underlying design principles. Section <ref> delves into the rigorous proof of preservation of convexity for ICLSTM and its integration into a convex LMPC framework. Section <ref> validates the performance and efficiency of our proposed framework against established baselines through a case study involving a continuous stirred tank reactor (CSTR). Finally, Section <ref> encapsulates our key findings and draws conclusions based on the results.

§ NONLINEAR LYAPUNOV-BASED MODEL PREDICTIVE CONTROL

§.§ Notation

In the following sections, we adopt the common notation style in deep learning community, where $x$ denotes the input unless otherwise stated, $y$ denotes the output, $g$ denotes the activation function, $F$ denotes a function, $L$ denotes the number of layers, $W$ denotes the network weights, $b$ denotes the bias term, and $z$ denotes the hidden layer output. Set subtraction is denoted by $``\backslash"$, that is, $A \backslash B := \{x \in \mathbb{R} \ | \ x \in A, x \notin B\}$. Composition is denoted by $``\circ"$. Euclidean norm of a vector is denoted by $|\cdot|$. Moreover, $f$ denotes the forget gate, $i$ denotes the input gate, $o$ denotes the output gate, $c$ denotes the cell state, and $h$ denotes the hidden state in an LSTM network.

§.§ Class of Systems

In this work, we consider the class of systems that can be represented by the following class of ordinary differential equation (ODE):

\begin{equation}\label{eq:classofsystems}
\dot{x} = F(x,u)
\end{equation}

where $x \in \mathbb{R}^n $ denotes the state vector, $u \in \mathbb{R}^m $ is the control action. $F : D \times U  \to \mathbb{R}^n $ is $\mathcal{C}^1 $, where $D \subset \mathbb{R}^n $ and $U \subset \mathbb{R}^m $ are compact and connected subsets that contain an open neighborhood of the origin, respectively.  Throughout this work, we also assume that $F(0,0) = 0$, so that the origin $(x, u) = (0,0)$ is an equilibrium point. Since first-principles models may not be available for complex real-world systems, our goal is to develop a neural network for the nonlinear system of Eq. <ref> and incorporate it into MPC while ensuring the desired computational efficiency and closed-loop stability.

§.§ LMPC Formulation

The LMPC scheme using a neural network model as the prediction model is given by the following optimization problem:

\begin{gather}
\label{eq9a}
\mathcal{L} = \min_{u \in S(\triangle)} \int_{t_k}^{t_{k+N}}J(\tilde{x}(t),u(t))dt \\
\label{eq9b}
s.t. \ \dot{\tilde{x}}(t) = F_{nn}(\tilde{x}(t),u(t)) \\
\label{eq9c}
u(t) \in U, \ \forall t \in [t_k,t_{k+N}) \\
\label{eq9d}
\tilde{x}(t_k) = x(t_k) \\
\label{eq9e}
V({\tilde{x}}(t)) < V({x}(t_k)), \ if \ x(t_k) \in \Omega_{\rho} \backslash \Omega_{\rho_{nn}}, \ \forall t \in [t_k,t_{k+N})
\end{gather}

where $\tilde{x}$ is the predicted state trajectory, $S(\triangle)$ is the set of piecewise constant functions with period $\triangle$, and $N$ is the number of sampling periods in the prediction horizon. The objective function $\mathcal{L}$ in Eq. (<ref>) incorporates a cost function $J$ in terms of the system states $x$ and the control actions $u$ (i.e., $|x|^2 + |u|^2$). The dynamic function $F_{nn}(\tilde{x}(t),u(t))$ in Eq. (<ref>) is parameterized as recurrent neural networks (e.g., plain RNN, ICRNN, ICLSTM, etc.). Eq. (<ref>) is the constraint function $U$ on feasible control actions. Eq. (<ref>) defines the initial condition $\tilde{x}(t_k)$  of Eq. (<ref>), which is the state measurement at $t = t_k$. Lastly, Eq. (<ref>) is the Lyapunov-based constraint $V$ that ensures closed-loop stability for the nonlinear system under LMPC by requiring that the value of $V(x)$ decreases over time, where $\Omega_{\rho}$ is the closed-loop stability region of the system, and $\Omega_{\rho_{nn}}$ is a small set around the origin where the state should ultimately be driven. The first element of the optimal input trajectory computed by the LMPC will be applied to the system over the sampling period and the LMPC is resolved again at the next sampling time. The detailed problem formulation, especially the design of the Lyapunov-based constraint, can be found in our previous works [Wu et al., 2019, Wu et al., 2019].

§ FAMILY OF INPUT CONVEX NEURAL NETWORKS

§.§ Input Convex Neural Networks

ICNNs represent a category of deep learning models where the output is designed to exhibit convexity with respect to the input. Currently, there exist two primary variants of input convex architectures: Input Convex Feedforward Neural Network (ICFNN) (i.e., Fig. <ref>) and Input Convex Recurrent Neural Network (ICRNN) (i.e., Fig. <ref>).

§.§.§ ICFNN

The output of each layer of ICFNN follows:

\begin{equation}
\label{eq_icfnn}
z_{t+1} = g_t(W_t^{(z)}z_t + W_t^{(x)}x + b_t)
\end{equation}

The output $z_{t+1}$ is input convex for single-step prediction if all weights $W_{1:L-1}^{(z)}$ are non-negative and all activation functions $g$ are convex and non-decreasing [Amos et al., 2017], while the output $z_{t+1}$ is input convex for multi-step ahead predictions if all weights $W_t^{(z)}$ and $W_t^{(y)}$ are non-negative and all activation functions $g$ are convex and non-decreasing [Bünning et al., 2021].

§.§.§ ICRNN

The output of ICRNN follows:

\begin{gather}
h_t = g_1(U\textbf{x}_t + Wh_{t-1} + D_2\textbf{x}_{t-1}) \\
y_t = g_2(Vh_t + D_1h_{t-1} + D_3\textbf{x}_t)
\end{gather}

The output $y_t$ is input convex if all weights $[U, W, V, D_1, D_2, D_3]$ are non-negative and all activation functions $g$ are convex and non-decreasing, where $\textbf{x}$ denotes the input $\big[\begin{smallmatrix}
x \\ -x
\end{smallmatrix}\big]$ [Chen et al., 2018].

However, the current ICNNs suffer from the vanishing gradient problem, which hinders their ability to effectively model complex data as deep neural networks. This limitation has motivated us to develop ICLSTM.

§.§ Introduction of Input Convex LSTM

In this work, we introduce a novel input convex architecture based on LSTM [Hochreiter and Schmidhuber, 1997], referred to as Input Convex LSTM and depicted in Fig. <ref>, which is designed to address the vanishing gradient problem of the current architectures while maintaining a similar level of model complexity, especially for time-series tasks.

Architecture of ICNNs

Specifically, the output of ICLSTM follows:

\begin{gather}
h_{lstm} = F(\textbf{x};\mathbf{W_t^{(lstm)}}) \\
z_t = g_1[W_t^{(z)}h_{lstm} + b_t^{(z)}] + \textbf{x} \\
y_t = g_2[W_t^{(y)}z_t + b_t^{(y)}]
\end{gather}

where $\textbf{x}$ denotes the input $\big[\begin{smallmatrix}
x \\ -x \end{smallmatrix}\big]$, $\theta = [\mathbf{W_t^{(lstm)}}, W_t^{(z)}, W_t^{(y)}]$ denotes the non-negative network weights, $g_1$ denotes the convex, non-negative, and non-decreasing activation functions, $g_2$ denotes the convex and non-decreasing activation function, $h_{lstm}$ denotes the output of the LSTM layer which is a function with respect to $\textbf{x}$ with a convex, non-negative, and non-decreasing initial hidden state $h_0$ and initial cell state $c_0$. If we unroll the dynamics with respect to time, we have $y_t = F(\mathbf{x_1},\mathbf{x_2},...,\mathbf{x_t};\theta)$, where $\theta = [\mathbf{W_t^{(lstm)}}, W_t^{(z)}, W_t^{(y)}]$.

Similar to ICRNN, we expand the input $x$ to $\big[\begin{smallmatrix} x \\ -x\end{smallmatrix}\big]$ by including its negation duplicates [Chen et al., 2018]. This concept initially aimed to facilitate network composition in dynamic system scenarios, but it yields additional advantages. In our experiment, we discovered that incorporating a non-negative weight constraint in an ICNN restricts its representability. The expanded input allows for a more accurate representation of dynamic systems as opposed to the original input. This input expansion can be regarded as a form of data augmentation, bolstering the model's robustness. Furthermore, including the negation of the input enhances gradient flow during training. By providing the network with both $x$ and $-x$, we introduce a greater number of symmetric data points, resulting in more consistent and well-balanced gradients throughout the training process.

While the proposed ICLSTM shares similarities with ICRNN and ICFNN,  it should be pointed out that there are major differences between the two, and that ICLSTM is not simply replacing RNN layers by LSTM layers. Specifically, unlike ICFNN and ICRNN, our ICLSTM enforces an additional non-negative constraint on the activation function and the initial hidden state and cell state, which is to ensure the convexity of the network. Furthermore, ICFNN and ICRNN leverage weighted direct “passthrough" layers to enhance representational capabilities, while ICLSTM adopts a parameter-free skip connection to mitigate the issues of vanishing gradients and enhance generalization [He et al., 2016].
Moreover,  a dense layer with the same dimension as the input is followed by every LSTM layer to maintain consistent dimensions between the input and output of the LSTM layer. This configuration facilitates the subsequent concatenation of the layer output with the input via the skip connection.  It is worth noting that the use of parameter-free skip connections and dense layers reduces network complexity. This simplification is particularly beneficial as it offsets the internal complexity inherent to the LSTM layer. In the following section, we will prove the convexity of the proposed ICLSTM network.

§ CLOSED-LOOP CONTROL USING ICLSTM-BASED MPC

The following lemmas are provided to prove Theorem <ref> on the convexity of a L-layer ICLSTM.

§.§ Conveixty of ICLSTM

If there are two or more convex functions, e.g., $f : \mathbb{R}^n \rightarrow \mathbb{R}$ and $g : \mathbb{R}^m \rightarrow \mathbb{R}$, by taking their weighted sum with non-negative coefficients, i.e., $\alpha f + \beta g$ where $\alpha \geq 0$ and $\beta \geq 0$, the result is still a convex function [Boyd and Vandenberghe, 2004].

Suppose $f : \mathbb{R}^n \rightarrow \mathbb{R}$, $A \in \mathbb{R}^{n \times m}$, and $b \in \mathbb{R}^n$. Define $g : \mathbb{R}^m \rightarrow \mathbb{R}$ by $g(x) = f(Ax+b)$, with $dom(g) = \{x~|~Ax+b \in dom(f)\}$. Then, if $f$ is convex, $g$ is convex [Boyd and Vandenberghe, 2004].

If there is a convex function $f : \mathbb{R}^n \rightarrow \mathbb{R}$ and by composing it with a monotone non-decreasing function $g : \mathbb{R} \rightarrow \mathbb{R}$, i.e., $h = g \circ f$, then $h$ is convex and non-decreasing [Boyd and Vandenberghe, 2004].

Let $f : \mathbb{R}^n \rightarrow \mathbb{R}$ and $g : \mathbb{R}^m \rightarrow \mathbb{R}$ be two convex monotone non-decreasing functions defined on convex sets. The product $h = f \times g$ is convex if and only if both $f$ and $g$ are non-negative on their respective domains.

Let $h = fg$, then $h' = gf'+fg'$ and $h'' = gf''+fg''+2f'g'$. A function is convex if and only if its second derivative is non-negative over its entire domain. Thus, a necessary and sufficient condition for $h$ to be convex is $h'' \geq 0$. Since $f$ and $g$ are convex, $f'' \geq 0$ and $g'' \geq 0$. Since $f$ and $g$ are monotone non-decreasing, $f'$ and $g'$ have the same sign. Assuming both $f$ and $g$ are non-negative, $f \geq 0$ and $g \geq 0$. By combining aforementioned inequalities, we show that $gf'' \geq 0$, $fg'' \geq 0$, $f'g' \geq 0$, and lastly $h'' \geq 0$. Hence, the necessary and sufficient condition for the product of two convex monotone non-decreasing function $h=fg$ remains convex is that both $f$ and $g$ have to be non-negative on their respective domain.

Next, we develop the following theorem to show the convexity of L-layer ICLSTM.

Consider L-layer ICLSTM as shown in Fig. <ref>, the output $y$ is a convex non-decreasing function of the elements of the input $\textbf{x}$, where $\textbf{x}$ denotes $[x, -x]$ and $x$ is in a convex feasible space $\mathcal{X}$ (i.e., $x \in \mathcal{X}$), if and only if all of the following three conditions are met: (1) All weights $W_i$ are non-negative; (2) All activation functions $g_i$ are convex, non-decreasing, and non-negative (e.g., ReLU), except for the activation function of the output layer $g^y$ which is convex and non-decreasing (e.g. Linear); (3) The initial hidden state $h_0$ and the initial cell state $c_0$ are convex, non-decreasing and non-negative (e.g., $h_0$ and $c_0$ are initialized as zero vectors).

The output of ICLSTM follows:

\begin{gather}
\label{eq1} f_1 = g_1^{f}[W_1^{f}h_0 + W_1^{f}\textbf{x} + b_1^{f}] \\
\label{eq2} i_1 = g_1^{i}[W_1^{i}h_0 + W_1^{i}x + b_1^{i}] \\
\label{eq3} c_1 = g_1^{c}[W_1^{c}h_0 + W_1^{c}x + b_1^{c}] \\
\label{eq4} o_1 = g_1^{o}[W_1^{o}h_0 + W_1^{o}x + b_1^{o}] \\
\label{eq5} c_{new} = f_1c_0 + i_1c_1 \\
\label{eq6} h_{new} = o_1g_1^{h}(c_{new}) = h_{lstm} \\
\label{eq7} z_1 = g_1^{d}[W_1^{d}h_{lstm} + b_1^{d}] + \textbf{x} \\
\label{eq8} y = g^{y}[W_1^{y}z_1 + b_1^{y}]
\end{gather}

Let us consider a 1-layer ICLSTM, for Eq. (<ref>), the elements of $W_1^{f}h_0$ and $W_1^{f}\textbf{x}$ are convex functions assuming that all elements of $W_1^{f}$ are non-negative and $h_0$ is convex in its inputs (i.e., Lemma <ref>). The bias term $b_1^{f}$ is a constant. Thus, the term $W_1^{f}h_0 + W_1^{f}\textbf{x} + b_1^{f}$ is convex in $\textbf{x}$ (i.e., Lemma <ref> and <ref>) and the term $g_1^{f}[W_1^{f}h_0 + W_1^{f}\textbf{x} + b_1^{f}]$ is convex, non-decreasing and non-negative in $\textbf{x}$ (i.e., Lemma <ref>). Without loss of generality, $i_1$, $c_1$ and $o_1$ are convex, non-decreasing and non-negative.

Since $f_1$, $c_0$, $i_1$, $c_1$ are non-negative non-decreasing convex functions, $c_{new}$ is convex (i.e., Lemma <ref> and <ref>), and thus $h_{new}$ is convex (i.e., Lemma <ref> and <ref>). Without loss of generality, $f_t$, $c_t$, $i_t$, $c_t$, $c_{new_t}$ and $h_{new_t}$ in the following $t$ iterations are convex. Up to this step, we have proved that the output of LSTM $h_{lstm}$ is input convex.  Subsequently, without loss of generality, $z_1$ is convex and thus $y$ is convex. Therefore, we have shown that the 1-layer ICLSTM output $y$ is a convex non-decreasing function of the elements of the input $\textbf{x}$.

In the $i^{th}$ layer, the input of the LSTM network $\textbf{x}$ is substituted with the output of previous layer $z_{i-1}$. Since $z_{i-1}$ is a convex, non-decreasing function with respect to the input $\textbf{x}$, without loss of generality, the L-layer ICLSTM output $y$ is a convex non-decreasing function of the elements of the input $\textbf{x}$.

§.§ Application of ICLSTM to a Finite-horizon Convex LMPC

Next, we utilize an ICLSTM to model the state transition dynamics, as expressed by $\dot{\tilde{x}}(t) = F_{nn}(\tilde{x}(t),u(t))$. This neural network is then embedded into a finite-horizon convex LMPC problem as designed in Eq. (<ref>). The primary objective of this integration is to determine the optimal sequence of actions, denoted as $u_t, u_{t+1}, \ldots, u_{t+N}$, for a predetermined  prediction horizon $N$. We first present a lemma to show the necessary and sufficient conditions for the designed LMPC in Eq. (<ref>) to be convex.

By embedding ICLSTM into the LMPC designed in Eq. (<ref>), the optimization problem is considered as a convex optimization problem if and only if both the objective function and constraints are convex.

Let us define the control Lyapunov function to be $V(x) = x^TPx$, where $x \in \mathbb{R}^{n \times 1}$ and $P \in \mathbb{R}^{n \times n}$. For $V(x)$ to exhibit convexity, it is necessary for its Hessian matrix, denoted as $H$, to be positive semidefinite. In this case, $H$ is equal to $2P$. Therefore, Eq. (<ref>), represented as $V(\tilde{x}(t)) - V(x(t)) < 0$, is convex according to Lemma <ref> without any loss of generality, by selecting a positive semidefinite matrix $P$.

Furthermore, the cost function of Eq. (<ref>) for a tracking MPC that is typically designed in a quadratic form, i.e., $|\tilde{x}(t)-0|^2+|u(t)-0|^2$ ($(0, 0)$ is the steady state)   is a convex function. Given that Eq. (<ref>) is parameterized as ICLSTM, and both Eq. (<ref>) and Eq. (<ref>) take the form of affine functions, the LMPC problem outlined in Eq. (<ref>) qualifies as a convex optimization problem, provided that the $P$ matrix is designed to be positive semidefinite.

Next, we develop the following theorem to show a convex LMPC with multi-step ahead prediction is still convex.

Consider a neural network-based convex LMPC, the problem remains input convex in the face of multi-step ahead prediction (i.e., when the prediction horizon $N$ is greater than 1), if the neural network embedded within the LMPC of Eq. (<ref>) is inherently input convex (e.g., ICLSTM).

The proof of Theorem <ref> is intuitive. Consider a $2$-step ahead prediction problem (i.e., N = 2) with a $n$-layer embedded ICLSTM $f_i(x_i,u_i)$, the final output is $y_2 = f_2(x_2 = f_1(x), u_2)$, where $\textbf{x}$ is the input. It is equivalent to a $1$-step ahead prediction problem with a $2n$-layer embedded ICLSTM, but with a new input $u_2$ concatenated at the $n^{th}$ layer. Without the loss of generality from Theorem <ref>, the $2$-step ahead prediction remains input convex. Hence, without loss of generality, a $m$-step ahead prediction problem with a $n$-layer embedded ICLSTM is equivalent to a $1$-step ahead prediction problem with a $mn$-layer embedded ICLSTM with new inputs $u_i$ concatenated at every $n^{th}$ layers, which is indeed input convex.

§ APPLICATION TO A CHEMICAL PROCESS

Finally, we apply LMPC to a well-mixed, nonisothermal continuous stirred tank reactor where an irreversible second-order exothermic reaction takes place (i.e., the reaction will transform a reactant A to a product B). The CSTR is equipped with a heating jacket that supplies/removes heat at a rate $Q$. The CSTR dynamic model is described by the following material and energy balance equations:

\begin{gather}
\frac{dC_A}{dt} = \frac{F}{V}(C_{A0}-C_A)-k_0e^{\frac{-E}{RT}}C_A^2 \\
\frac{dT}{dt} = \frac{F}{V}(T_0-T)+\frac{-\triangle H}{\rho_LC_p}k_0e^{\frac{-E}{RT}}C_A^2+\frac{Q}{\rho_LC_pV}
\end{gather}

where $C_A$ is the concentration of reactant A, $V$ is the volume of the reacting liquid, $T$ is the temperature, $Q$ is the heat input rate, $C_{A0}$ is the inlet concentration of reactant A, $T_0$ is the inlet temperature, $F$ is the volumetric flow rate, $\rho_L$ is the constant density of the reacting liquid, $C_p$ is the hear capacity, $\triangle H$ is the enthalpy of reaction, $k_0$ is the pre-exponential constant, $E$ is the activation energy, and $R$ is the ideal gas constant.

The manipulated inputs in this system are represented by $\triangle C_{A0} = C_{A0} - C_{A0_s}$ and $\triangle Q = Q - Q_s$, which correspond to the inlet concentration of reactant A and the heat input rate, respectively. The states of the closed-loop system can be described as $x^T = [C_A - C_{As}, T - T_s]$, while the control actions are denoted as $u^T = [\triangle C_{A0}, \triangle Q]$, such that the equilibrium point of the system is located at the origin of the state-space. The main control objective is to operate the CSTR at the unstable equilibrium point $(C_{As}, T_s)$ by manipulating $\triangle C_{A0}$ and $\triangle Q$, using the LMPC in Eq. (<ref>) with neural networks, and finally reach the steady-state.

We first perform open-loop simulations in the closed-loop stability region $\Omega_\rho$ for the CSTR of Eq. (<ref>) to train the neural networks. The LMPC problem is then solved using PyIpopt, which is the python version of IPOPT [Wächter and Biegler, 2006], with an integration time step of $h_c = 10^{-4}~hr$ and the sampling period $\triangle = 0.005~hr$.
The control Lyapunov function $V(x) = x^TPx$ is designed with the following positive definite $P$ matrix:

\begin{equation}
P = \begin{bmatrix}
1060 & 22 \\
22 & 0.52
\end{bmatrix}
\label{eq20}
\end{equation}

which ensures the convexity of the LMPC. The detailed system description and parameter values can be found in our previous works [Wu et al., 2019, Wu et al., 2019].

§.§ Control Performance

We constructed and trained two-layer neural networks, which are meticulously configured with a batch size of 256, the Adam optimizer, and the Mean Square Error (MSE) loss function. The model designs are shown in Table <ref>. The primary purpose of neural networks is to capture and encapsulate the system dynamics, subsequently integrating into the LMPC framework shown in Eq. (<ref>).

As the efficacy of LMPC hinges on two critical factors, namely the final state of convergence and the time taken to reach convergence, relying solely on the training and validation MSE of neural networks may not offer a balanced assessment of their performance. To rectify this, we introduce a new performance metric specifically tailored for the LMPC task. This metric places a specific emphasis on the temporal aspect by evaluating the time taken for the neural network-based LMPC to achieve convergence. For the CSTR example, we define the small region of convergence to be $|C_{A_i}| < 0.1$ and $|T_i| < 3$, where $C_{A_i} = C_A - C_{As}$ and $T_i = T - T_s$, and the equation of the stability region is defined as $1060x^2 + 44xy + 0.52y^2 - 372 = 0$, which is an ellipse in state space.

Model     Activation     Weight constraint     No. of layers     No. of neurons     Validation MSE     No. of parameters

Plain RNN     Tanh     No constraints     2     64     2.7492e-05     12,802

Plain LSTM     Tanh     No constraints     2     64     1.7523e-06     50,818

ICFNN     ReLU     Non-negative     2     64     0.1531     13,076

ICRNN     ELU     Non-negative     2     64     4.7460e-04     38,530

ICLSTM     ReLU     Non-negative     2     64     2.4660e-06     38,434

Hyperparameter design of models

The validation MSE presented in Table <ref> reveals a suboptimal performance of ICFNN for LMPC-based CSTR which can be characterized as a time-series forecasting task. ICFNN struggles to effectively capture the system dynamics, which will eventually lead to the divergence of LMPC. Consequently, we omit ICFNN from the subsequent comparative analysis.

§.§.§ Vanishing gradient problem on ICRNN

Training losses on ICRNN and ICLSTM w.r.t. number of layers

The training losses, measured in MSE, for both ICRNN and ICLSTM are visualized in Fig. <ref>. As the depth of the neural network increases, the training loss of ICRNN also increases. This behavior suggests the presence of a vanishing gradient problem, limiting the deep ICRNN's capacity to effectively capture and generalize system dynamics.

In contrast, the training loss for the proposed ICLSTM remains consistently low and close to zero as the number of layers increases. This outcome highlights the ICLSTM's ability to learn and adapt to complex dynamics more effectively when compared to ICRNN, owing to its inherent capacity as a deep neural network.

§.§.§ Control performance of ICLSTM

Our experimentation focuses on evaluating the control performance of the CSTR, as defined by Eq. (<ref>), by embedding ICLSTM to LMPC designed in Eq. (<ref>). To ensure consistency and reliability in our evaluation, we deliberately restrict the selection of random initial conditions, denoted as $(C_{Ai}, T_i)$, to be distant away from the system's equilibrium point. This measure aims to minimize the influence of computational randomness. When the initial conditions are situated near the equilibrium point, the system converges rapidly within a few iterations, which makes the purpose of comparison in convergence runtime insignificant.

2|c|Plain RNN     2|c|Plain LSTM     2|c|ICRNN     ICLSTM

$\mathbf{[C_{A_i}, T_i]}$    Time (s)     % Decrease     Time (s)     % Decrease     Time (s)     % Decrease     Time (s)     Better?

$[-1.5, 70]$     2201.86     53.730%     1708.84     40.381%     1204.34     15.407%     1018.79     $\surd$

$[1.5, -70]$     2554.50     39.765%     1938.92     20.641%     1733.70     11.248%     1538.70     $\surd$

$[-1.3, 60]$     1412.78     54.058%     1076.55     39.709%     1114.71     41.773%     649.06     $\surd$

$[1.3, -60]$     2194.86     31.073%     1552.72     2.567%     1403.33     -7.805%     1512.86     $\times$

$[-1, 55]$     1259.34     25.561%     1450.84     35.386%     1065.07     11.983%     937.44     $\surd$

$[1, -55]$     1767.98     44.210%     1106.84     10.886%     1037.69     4.948%     986.35     $\surd$

$[-1.25, 50]$     1312.55     24.292%     1327.58     25.149%     1246.42     20.275%     993.71     $\surd$

$[1.25,-50]$     1974.38     64.540%     1854.80     35.306%     1353.39     11.338%     1199.94     $\surd$

$[-0.75, 40]$     1787.51     60.532%     1174.78     39.947%     1329.34     46.929%     705.49     $\surd$

$[0.75, -40]$     1923.91     69.501%     1584.45     62.966%     1090.99     46.216%     586.78     $\surd$

Average     1839.0     46.7%     1477.6     31.3%     1258.0     20.2%     1012.9     $\surd$

Convergence runtime of neural network-based LMPC and its respective percentage decrease with respect to ICLSTM-based LMPC

In this experiment, all trials successfully achieved convergence. Table <ref> presents compelling results, showcasing that ICLSTM-based LMPC yields an improvement in convergence runtime. Specifically, it attains an average percentage decrease of 46.7%, 31.3%, and 20.2% compared to plain RNN, plain LSTM, and ICRNN, respectively.

§ CONCLUSION

In this study, we introduced a novel neural network architecture engineered to ensure convexity of the output with respect to the input, specifically tailored for finite-horizon convex LMPC. Notably, our framework excels in terms of control performance, showcasing improvements in convergence runtime. Furthermore, our framework enables deep neural network training by successfully addressing the vanishing gradient problem, a hurdle that had previously hampered the effectiveness of prior ICNNs. Through the simulation study of a dynamic CSTR system, we demonstrated the efficacy and efficiency of our proposed framework. This work serves as a pivotal bridge between ICNNs and control applications within the realm of chemical engineering.

[Afram et al., 2017]
Abdul Afram, Farrokh Janabi-Sharifi, Alan S Fung, and Kaamran Raahemifar.
Artificial Neural Network (ANN) Based Model Predictive Control (MPC) and Optimization of HVAC Systems: A State of the Art Review and Case Study of a Residential HVAC System.
Energy and Buildings, 141:0 96–113, 2017.

[Alvarez-Melis et al., 2021]
David Alvarez-Melis, Yair Schiff, and Youssef Mroueh.
Optimizing Functionals on the Space of Probabilities with Input Convex Neural Networks.
arXiv preprint arXiv:2106.00774, 2021.

[Amos et al., 2017]
Brandon Amos, Lei Xu, and J Zico Kolter.
Input Convex Neural Networks.
In International Conference on Machine Learning, pages 146–155. PMLR, 2017.

[Bao et al., 2017]
Xuefeng Bao, Ziyue Sun, and Nitin Sharma.
A Recurrent Neural Network Based MPC for a Hybrid Neuroprosthesis System.
In 2017 IEEE 56th Annual Conference on Decision and Control (CDC), pages 4715–4720. IEEE, 2017.

[Boyd and Vandenberghe, 2004]
Stephen P Boyd and Lieven Vandenberghe.
Convex Optimization.
Cambridge University Press, 2004.

[Bünning et al., 2021]
Felix Bünning, Adrian Schalbetter, Ahmed Aboudonia, Mathias Hudoba de Badyn, Philipp Heer, and John Lygeros.
Input Convex Neural Networks for Building MPC.
In Learning for Dynamics and Control, pages 251–262. PMLR, 2021.

[Chen et al., 2018]
Yize Chen, Yuanyuan Shi, and Baosen Zhang.
Optimal Control via Neural Networks: A Convex Approach.
arXiv preprint arXiv:1805.11835, 2018.

[Chen et al., 2020]
Yize Chen, Yuanyuan Shi, and Baosen Zhang.
Data-Driven Optimal Voltage Regulation using Input Convex Neural Networks.
Electric Power Systems Research, 189:0 106741, 2020a.

[Chen et al., 2020]
Yize Chen, Yuanyuan Shi, and Baosen Zhang.
Input Convex Neural Networks for Optimal Voltage Regulation.
arXiv preprint arXiv:2002.08684, 2020b.

[Ellis and Chinde, 2020]
Matthew J Ellis and Venkatesh Chinde.
An Encoder-Decoder LSTM-Based EMPC Framework Applied to a Building HVAC System.
Chemical Engineering Research and Design, 160:0 508–520, 2020.

[He et al., 2016]
Kaiming He, Xiangyu Zhang, Shaoqing Ren, and Jian Sun.
Deep Residual Learning for Image Recognition.
In Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition, pages 770–778, 2016.

[Hochreiter and Schmidhuber, 1997]
Sepp Hochreiter and Jürgen Schmidhuber.
Long Short-Term Memory.
Neural Computation, 90 (8):0 1735–1780, 1997.

[Lanzetti et al., 2019]
Nicolas Lanzetti, Ying Zhao Lian, Andrea Cortinovis, Luis Dominguez, Mehmet Mercangöz, and Colin Jones.
Recurrent Neural Network Based MPC for Process Industries.
In 2019 18th European Control Conference (ECC), pages 1005–1010. IEEE, 2019.

[Makkuva et al., 2020]
Ashok Makkuva, Amirhossein Taghvaei, Sewoong Oh, and Jason Lee.
Optimal Transport Mapping via Input Convex Neural Networks.
In International Conference on Machine Learning, pages 6672–6681. PMLR, 2020.

[Nubert et al., 2020]
Julian Nubert, Johannes Köhler, Vincent Berenz, Frank Allgöwer, and Sebastian Trimpe.
Safe and Fast Tracking on a Robot Manipulator: Robust MPC and Neural Network Control.
IEEE Robotics and Automation Letters, 50 (2):0 3050–3057, 2020.

[Sitapure and Kwon, 2022]
Niranjan Sitapure and Joseph Sang-Il Kwon.
Neural Network-Based Model Predictive Control for Thin-film Chemical Deposition of Quantum Dots using Data from a Multiscale Simulation.
Chemical Engineering Research and Design, 183:0 595–607, 2022.

[Wächter and Biegler, 2006]
Andreas Wächter and Lorenz T Biegler.
On the implementation of an interior-point filter line-search algorithm for large-scale nonlinear programming.
Mathematical programming, 106:0 25–57, 2006.

[Wu et al., 2019]
Zhe Wu, Anh Tran, David Rincon, and Panagiotis D Christofides.
Machine-Learning-Based Predictive Control of Nonlinear Processes. Part i: Theory.
AIChE Journal, 650 (11):0 e16729, 2019a.

[Wu et al., 2019]
Zhe Wu, Anh Tran, David Rincon, and Panagiotis D Christofides.
Machine-Learning-Based Predictive Control of Nonlinear Processes. Part ii: Computational Implementation.
AIChE Journal, 650 (11):0 e16734, 2019b.

[Yang and Bequette, 2021]
Shu Yang and B Wayne Bequette.
Optimization-Based Control using Input Convex Neural Networks.
Computers & Chemical Engineering, 144:0 107143, 2021.

[Zhang et al., 2021]
Ling Zhang, Yize Chen, and Baosen Zhang.
A Convex Neural Network Solver for DCOPF with Generalization Guarantees.
IEEE Transactions on Control of Network Systems, 90 (2):0 719–730, 2021.

[Zheng et al., 2022]
Yingzhe Zheng, Xiaonan Wang, and Zhe Wu.
Machine Learning Modeling and Predictive Control of the Batch Crystallization Process.
Industrial & Engineering Chemistry Research, 610 (16):0 5578–5592, 2022a.

[Zheng et al., 2022]
Yingzhe Zheng, Tianyi Zhao, Xiaonan Wang, and Zhe Wu.
Online Learning-Based Predictive Control of Crystallization Processes under Batch-to-Batch Parametric Drift.
AIChE Journal, 680 (11):0 e17815, 2022b.