Flatness-Based Model Predictive Constrained Optimal Control for Chemical Reactor

Controlling a system with nonlinear behavior to achieve optimal performance is the ultimate challenge in engineering and practical sciences. Those systems' unpredictable performance and inherent complexity usually rebel against traditional control methods. Understanding and applying optimal control in nonlinear systems is vital for maintaining the stringent limits of system performance and dependability. Solving the Hamilton-Jacobi-Bellmen (HJB) equation has many difficulties, especially with complicated, multi-variable systems. These challenges motivate the pursuit of additional solutions to acquire optimal control.

In the field of process control, a continuous stirred tank reactor CSTR, described extensively in a study [1], presents an operational challenge due to its nonlinearity, output multiplicity, and oscillation [2]. In addition, CSTR has a versatile dynamic representation that depends mainly on the operating conditions and unstable equilibrium points [3]. The effects of disturbances and interaction between process variables on CSTR results in changes in the nominal operation region. Classical linear control is incompetent in regulating and trajectory-tracking problems in process systems, particularly in CSTR [4]. Consequently, advanced control techniques are popular in addressing the challenges of CSTR complex behaviors. At the top of these advanced techniques are predictive control [5], nonlinear control [6], robust control [7, 8], or a blend of these techniques as in the study [9]. The global transformation technique is essential for altering the original nonlinear CSTR model, using a new set of variables and state feedback, to an equivalent linear model [10]. Feedback linearization in the study [11] with an extended state observer was applied to the nonlinear CSTR. A technique for tracking a desired temperature utilizing the differential flat CSTR and feedforward control was introduced [12].

Differential flatness is an important feature suitable for controlling and planning trajectories for nonlinear systems. A system is differentially flat if an algebraic expression exists from its state and inputs by defining a set of flat outputs and their time derivatives [13]. This property facilitates the transformation of a possibly nonlinear control problem into a linear one using feedforward linearization, simplifying the design of control strategies [14]. Recent studies emphasize the successful implementation of flatness control in various situations where other control strategies failed or failed to give the desired performance [15]. Underactuated mechanical systems were considered, and simple flatness conditions were derived for those types of systems. In a study [16], the control of gas-turbine electric power generation units was performed using differential flatness theory. The flexible mobile manipulator trajectory tracking problem was discussed [17], and the flatness property was combined with active disturbance rejection control for the complicated control problem. Other successful stories of combining deferential flatness with other control techniques that suit special performance requirements, such as optimality, can be found extensively in the literature. An optimal control action can be obtained by combining differential flatness and linear quadratic regulator LQR as in studies [18, 19], resulting in superior performance compared to using LQR. The main restriction on using LQR is its inability to handle constraints imposed on the system's inputs.

For systems with critical constraints, techniques like Model Predictive Control (MPC) are often preferred, as they explicitly incorporate these constraints into the control design. MPC is an important control technology with significant influence in practice and an expected 85% future impact in many categories of industrial applications. This highly impactful percentage placed model predictive control in an advanced position ahead of many other control technologies, such as robust, adaptive, and nonlinear control [20]. MPC has the highest current impact in different fields, such as industry [21], aerospace [22], and automotive [23].

The analysis of MPC and its design was introduced in a framework that provides a diverse toolkit to apply and future development of MPC formulation [24]. The receding horizon control, or MPC, is an advanced control and optimization-based method. Its working principle can be explained as solving online, optimal control problems in an open-loop manner repeatedly within a finite horizon. Feedback is implicitly generated by implementing the first part of the optimized manipulated variable trajectory; this procedure is repeated at each time step. MPC comes with many advantages: (i) constraints consideration of state and control variables; (ii) applicability on both SISO or MIMO systems; (iii) utilization of general performance criterion for optimization.

Meanwhile, the MPC property of repeatedly solving an optimization problem online requires high computational effort. Although the number of transistors in a microprocessor has increased significantly, which lifts some of the computational burdens, the clock speed is limited due to heat dissipation in the transistors [25]. Another factor that increases the computational load is related to advantage (i), as imposing the constraints will affect the optimization problem that needs to be solved online. For systems with nonlinear dynamics in their intended operation region, using NMPC is suitable for solving the nonlinear control problem [26]. Granting that NMPC is an influential tool for controlling nonlinear systems, it is not always the best choice due to the requirement of solving nonlinear programming problems, which adds computational complexity in real time.

To simplify this problem, many solutions have been introduced to design a linear MPC with quadratic programming (QP) to solve the online optimization problem. Researchers adopted many approaches, including full or partial state feedback linearization FBL, to linearize the nonlinear system and use it with MPC in applications [27, 28]. Implementing a nonlinear state observer is mandatory to gain knowledge of the system state using this approach. Another issue concerning this method is the introduction of additional feedback, which contradicts the MPC structure. In addition, constraints on the inputs and states represent an added challenge to any control system, which can be efficiently handled by incorporating these limitations in linear MPC.

This paper proposes a joint control approach that consists of a flat CSTR system and linear MPC to solve the constrained optimization problem. This innovative formulation allows us to restate the MPC problem and turn it into a linear one, simplifying the control problem. The approach is designed to improve performance and constraint handling through a flatness-based model predictive control FMPC strategy. This control policy effectively maintains system performance and dependability and is relatively inexpensive from a computational perspective. It considers the importance of respecting the constraints on the input and output of the plant, providing a practical and cost-effective solution.

Integrating feedback MPC with Feedforward DF leverages the strengths of both approaches to create a more robust and efficient control strategy for nonlinear systems. The benefits of joining these control techniques can be listed as:

• The DF can reformulate the nonlinear control problem in a globally linearized flat space, allowing MPC to optimize over a simplified transformed system, which leads to faster convergence.
• FMPC provides an improved control performance due to the initial optimal trajectory generated by DF and refined by MPC in real-time, which balances optimality and feasibility.
• The integration of DF and MPC allows the controller to follow the desired trajectory and ensures compliance with system constraints. This enhances the robustness and uncertainty since the MPC adjusts accordingly in real time and adheres to constraints.

The rest of this article is structured as follows: Section 2 describes the CSTR dynamic model. Section 3 presents background knowledge on flatness and its feedforward linearization and a fundamental background on model predictive control. Section 4 presents the methodology. Section 5 presents the Simulation results. Section 6 concludes this article.

This section presents flatness with associated feedforward linearization and the fundamental design of MPC. Both techniques are used to design control algorithms separately and apply them to the nonlinear plants of CSTR. The procedure of combining the two control algorithms to form an FMPC strategy will be introduced as follows:

3.1 Differentially flat systems

Fliess et al. [30] introduced the concept known as the flatness of nonlinear systems. The important characteristic of flat systems is the existence of a finite set of variables called flat output or linearizing output. To give a brief background about differentially flat systems and their essential feature, we consider a nonlinear system presented by:

$\dot{x}(\mathrm{t})=f(x(t), u(t))$       (4)

$y(t)=h(x)$     (5)

with $x(t) \in R^n, u(t) \in R^m, y(t) \in R^m$ and $f, h$ are smooth functions in their domain. System (1) is said to be a differentially flat system if a flat output variable $\zeta$ exists and meets the following conditions:

Flat output $\zeta$ can be expressed as a function of the state:

$\zeta=\rho(x)$      (6)

The system state and input may be parametrized by the flat output and a finite number of its components and derivatives as a differential function.

$x=\psi_1\left(\zeta, \dot{\zeta}, \ldots, \zeta^{(n-1)}\right)$       (7)

$u=\psi_2\left(\zeta, \dot{\zeta}, \ldots, \zeta^{(n)}\right)$      (8)

where, $\psi_1$ and $\psi_2$ are smooth functions and  is the flat output required to describe the system. If the system complies with these requirements, the output of the actual system y is a function of the flat output $\zeta$ and finite number of its time derivatives:

$y=\psi_3\left(\zeta, \dot{\zeta}, \ldots, \zeta^{(n-\beta)}\right)$      (9)

where, $\beta$ is the relative degree. Using the flatness property represented by feedforward control law Eq. (8) and reference trajectory Eq. (7) in terms of the flat output and its derivative would provide the ability to obtain the variable of the original system's time function without integration. In general, the differentially flat system can be represented using Brunovsky or flat states:

$\boldsymbol{z}:=\left[\zeta_1, \dot{\zeta}_1, \ldots, \zeta_1^{\left(\beta_1-1\right)}, \ldots, \zeta_m, \ldots, \zeta_1^{\left(\beta_m-1\right)}\right]$         (10)

$\zeta_i^{\left(\beta_i\right)}=\alpha_i\left(\zeta_1, \dot{\zeta}_1, \ldots, \zeta_1^{(\beta-1)}, u, \dot{u}, \ldots, u^{\left(\beta_i\right)}\right):=v_i$        (11)

where, the highest derivative $\beta_i$of the flat output $\zeta$ considered as the new flat input $v$. For multi-input flat system:

$\mathbf{v}:=\left[v_1, v_1, \ldots, v_m\right]^T$      (12)

Note that, in Eq. (11), $\beta_i$ is related to the relative degree of the system, while the number of flat inputs $v_i$ depends on the number of selected flat outputs $\zeta_i$ that represent the nonlinear system. Standard form can be obtained from (10) and (12) as:

$\dot{\mathbf{z}}=A \boldsymbol{z}+B v$       (13)

$\mathbf{v}=\psi_4^{-1}\left(z, u, \dot{u}, \ldots, u^{(\beta)}\right)$        (14)

where, (13) can be termed as a linear flat model. The method of feedforward linearization allows trajectory generators or controllers, such as MPC in our proposal, to only consider the linear flat model with the condition of matching initials for both the desired output and the flat output. To this end, we can conclude the section with the following theorem:

Theorem 1: (14) Considering a desired trajectory in the flat output $\zeta_d$ and a desired flat input $v_d$, with initial condition $z(0)=z_d(0)$, applying the nominal feedforward control law

$u=\psi_2\left(z_d, v_d\right)$       (15)

to the flat system, Eq. (4), will result in an equivalent linear system in Brunovsky canonical form by change of coordinates.

Using Theorem 1 allows controllers, such as MPC in our proposed method, to consider a linearized flat model. The output of MPC can then be fed to the nonlinear system through the inverse transformation Eq. (15) to correct the nonlinear relationship in Eq. (4). Using the desired flat state and input, instead of the actual measured flat state, represents the main difference between feedforward and feedback linearization.

3.2 Model predictive control

Predictive control determines an optimal future sequence of the manipulated variable. This procedure aims to minimize a cost function concerning the controlled variable future behavior on a finite time interval called the prediction horizon. The method of receding horizon, which can be described as optimizing at each time instant and only executing the first value of the obtained manipulated variable sequence, is utilized. The future behavior of the systems can be predicted based on actual measurements of a system model. This model is in state space representation, and the MPC problem is formulated in discrete time. Offset-free tracking, in the presence of constant disturbances or mismatch in the plant model, can be achieved by the augmented state space model:

$\left[\begin{array}{c}x(k+1) \\ d_c(k+1)\end{array}\right]=\left[\begin{array}{cc}A & B_{a c} \\ 0 & I\end{array}\right]\left[\begin{array}{c}x(k) \\ d_c(k)\end{array}\right]+\left[\begin{array}{c}B \\ 0\end{array}\right] u(k)$                   (16)

$y(k)=\left[\begin{array}{ll}C & C_{d c}\end{array}\right]\left[x(k) \quad d_c(k)\right]^T$        (17)

where, $x \in \mathbb{R}^n$ is the state vector, $u \in \mathbb{R}$ is the input variable, $d_c \in \mathbb{R}$ the disturbance and $y \in \mathbb{R}$ is the output. $A, B, C$ are matrices with appropriate scales and are known. The matrices of the disturbances, namely; $B_{d c}$ and $C_{d c}$ are designed based on the knowledge of disturbance. Note that the system is assumed to be stable, controllable, and observable.

The nonlinear plant is assumed to have the following input constraints on $u$:

$u_{\min } \leq u \leq u_{\max }, \Delta u_{\min } \leq \Delta u \leq \Delta u_{\max }$        (18)

A cost function $J(k)$ is need to be minimized, formulated as quadratic form, must consider the future deviation of the controlled variable from its desired steady state:

$\begin{array}{r}J(k)=\sum_{i=1}^{H_p}\left(\hat{y}(k+i)-C x_{s s}(k)\right)^T Q(\hat{y}(k+i) \left.-C x_{s s}(k)\right)+(\widehat{u}(k+i-1) \left.-u_{s s}(k)\right)^T R\left(\widehat{u}-u_{s s}(k)\right)\end{array}$                (19)

where, $\hat{y}(+i)$ is the predicted output at the time $k+i$, evaluated by measurements at the time $k$, and $\widehat{u}(k+i)$ is the predicted input. $H_p$ is the prediction horizon length and $Q, R$ are design parameters, $x_{s s}, u_{s s}$ the steady-state values of the state and the input can be obtained from Eq. (16). A linear disturbance observer is necessary to estimate the constant disturbance value. This observation closes the control loop similarly to the internal model control structure. The difference between actual plant and internal model output and the measure for disturbance and plant model mismatch is used to shift the value of the controlled variable. If the constraint on the input variable is applied, then the inequality constraint is a condition for the minimization of the cost function:

$u_{\min }(k+i-1) \leq u(k+i-1) \leq u_{\max }(k+i-1)$        (20)

where, $u_{\min }, u_{\max }$ represent the input variable limitations. Setting the disturbance to zero, the problem of minimizing the cost function subjected to constraints can be reformulated as:

$\min _u J(k)=\frac{1}{2} \mathcal{U}(k)^T \mathcal{H} U(k)-\mathcal{G}^T U(k)$             (21)

Subject to: $\Omega \mathcal{U}(k) \leq \omega(k)$      (22)

where, $U(k) \in \mathbb{R}^{H_p \times 1}$ has the input sequence within the prediction horizon. $\mathcal{H}, \mathcal{G}$ are known matrices that are appropriately scaled, $\Omega, \omega(k)$, represent the constraints. The optimization problem, which is composed of Eq. (21) and Eq. (22) and constraints Eq. (20), is a standard quadratic programming (QP) problem.

This section presents the simulation results of the FMPC, FFC, and LMPC controllers applied to the nonlinear CSTR system. To evaluate the performance of tracking a reference temperature trajectory using our proposed method and other control strategies. The second-order CSTR was employed with an unstable operating point and system parameters provided in Table 1. Eqs. (1)-(3) are used to simulate the dynamics of the CSTR. The manipulated input $u$ was constrained with bound constraints as $280 \leq T_c \leq 380$, no output constraints were imposed in this study. The comparison was made among three controllers, namely, linear model Predictive control LMPC, flatness feedforward control FFC, and compared against flatness-based model predictive control FMPC.

The description of each control method used in the simulation is as follows:

• The LMPC schematic diagram is shown in Figure 4, using an augmented linear model as described in Eqs. (16) and (17). It was implemented with tuning parameters, $H_p$ = 10, $H_c$ = 2, MV weight rate = 0.1, OV weight = 1. Three inputs have been selected in the design of the LMPC internal model, $T_C$ = MV, $T_f$ = MD and $C_{A f}$ = UD. The augmented system was trimmed to obtain a linear model and used reactor temperature $T$ = MO and $C_A$ = UO. The closed-loop performance in terms of the prediction model and the feedback connection of the nonlinear plant with the controller was tuned to be stable.

4.png

Figure 4. Linear model predictive control schematic

5.png

Figure 5. Flatness feedforward control schematic

FFC was designed to track the trajectory of nonlinear CSTR, as shown in Figure 5. Using the flatness property, trajectory tracking is reduced to a problem of a linear system required to track a given trajectory. Imposing an asymptotic stable behavior to the flat output tracking error, designed based on linear time-invariant dynamics, is the most straightforward approach to deal with nonlinear system tracking problems. By giving a desired trajectory $\zeta^*(t)$ which is required to be tracked by the flat output $\zeta$ the error dynamics can be defined as:

• $e(t)=\zeta-\zeta^*(t)$. The following linear feedback controller was used:

$v=\frac{d^\alpha}{d t^\alpha} \zeta^*(t)-\sum_{j=0}^{\alpha-1} \gamma_j \frac{d^{\alpha-j-1}}{d t^{\alpha-j-1}}\left(\zeta-\zeta^*(t)\right)$          (36)

The error evolved to stability by apparated choice of the coefficients of the linear feedback controller $\gamma_j$. The characteristic polynomial of the closed-loop system is Hurwitz, and the error signal was exponentially asymptotically stable to zero.

• To implement FMPC optimal control strategy, we started by utilizing the nonlinear CSTR differential flatness property; we defined the flat output state $z$ Eqs. (24), (25) and input $v$ Eq. (27), then we rewrote the nonlinear model as a discrete linear flat model as in Eq. (28). An MPC formulation was used subjected to a discretized linear flat model, Sampling time = 0.05, prediction horizon =10, and control horizon = 2. The schematic of the proposed method is presented in Figure 3. FMPC, described in detail earlier, exhibited high performance in tracking the setpoint with less control efforts when compared to the other control techniques.

Three simulation scenarios were investigated to demonstrate the effectiveness of the proposed control strategy in testing the performance of FMPC. The performance of the three controllers is compared for the setpoint change of $\pm$25 K from its nominal value, Figure 6 and Figure 7 are for the positive setpoint change.

6.png

Figure 6. Positive setpoint change comparison for reactor temperature

7.png

Figure 7. Positive setpoint change comparison for cooling temperature

As shown in Figure 7, the three controllers had to follow a positive setpoint change of 25 K above their nominal value of 350 K. The same bound on the control input was applied to this case, and only the output of FMPC showed no overshoot in its response. Table 2 presents the numerical values of the transient response analyses obtained using the stepinfo command in MATLAB.

Table 2. Transient response characteristics for positive setpoint change

Table 3 compares the parameters of the performance indicator and quantifies the control performance. The ISE index rates the efficiency of the three controllers in holding on to the output variable near the reference. ITAE and ITSE indicate the oscillation in the error signal. RMSE represents the average magnitude of prediction errors.

Table 3. Performance indicators for positive setpoint change

The result showed a 46.45% improvement in FMPC compared to LMPC (RMSE).

Figures 8 and 9 show negative setpoint change simulation results, with the reference setpoint temperature dropping to 325 K from 350 K.

8.png

Figure 8. Negative setpoint change comparison for reactor temperature

9.png

Figure 9. Negative setpoint change comparison for cooling temperature

The FMPC has shown a response capable of regulating the temperature with less overshoot. FFC and FMPC have comparable transient response results, while LMPC has a sluggish response due to its inaccurate linearized model and poor constraints handling.

The performance of the three controllers shown in Figure 10 and Figure 11 are the results of closed-loop response for feed-temperature $T_f$ disturbance. The input feed temperature was increased by 45 K while keeping the other nominal values as in Table 1, and the input constraints were still applied. LMPC response was prolonged and showed significant control moves, while FFC and FMPC had comparable results with the advantage of faster response to FMPC.

10.png

Figure 10. Feed temperature disturbance comparison reactor temperature

11.png

Figure 11. Feed temperature disturbance comparison reactor temperature

Table 4. A comparison of the elapsed computational time

The three controllers' simulations used a sampling time of 50 ms to measure the computation time of the controlled CSTR system. Table 4 includes the Average Computation Time (A_CT) and Maximum Computation Time (M_CT) of the LMPC, FFC, and FMPC controllers. The results show an enhanced average computation time of the FMPC compared to the other controllers. The real-time feasibility analysis showed that the average and worst-case computation times are well below the control loop sampling time with the advancement of FMPC, allowing implementation feasibility for more complicated systems. The results of Table 4 have been obtained based on a built-in timing function to capture the computation time of each of the applied control loops.

To conclude this section, FMPC demonstrates the ability to regulate the reactor's temperature setpoint while considering the constraints on the manipulated variable, and it shows good disturbance rejection. In addition, FMPC delivers the best results for both positive and negative setpoint changes and disturbance rejection, while LMPC produces the least favorable response.