Model Predictive Control with a Relaxed Cost Function for Constrained Linear Systems

-e Model Predictive Control technique is widely used for optimizing the performance of constrained multi-input multi-output processes. However, due to its mathematical complexity and heavy computation effort, it is mainly suitable in processes with slow dynamics. Based on the Exact Penalization -eorem, this paper presents a discrete-time state-space Model Predictive Control strategy with a relaxed performance index, where the constraints are implicitly defined in the weighting matrices, computed at each sampling time.-e performance validation for theModel Predictive Control strategy with the proposed relaxed cost function uses the simulation of a tape transport system and a jet transport aircraft during cruise flight. Without affecting the tracking performance, numerical results show that the execution time is notably decreased compared with two well-known discrete-time state-space Model Predictive Control strategies. -is makes the proposed Model Predictive Control mainly suitable for constrained multivariable processes with fast dynamics.


Introduction
Model Predictive Control (MPC) for linear systems is now a well-established discipline providing stability, feasibility, and robustness [1][2][3][4][5][6]. Due to its inherent ability to take into account constraints and deal with multi-input multi-output variables [7][8][9][10], it has been applied in a wide range of applications, including chemical processes, industrial systems, energy, health, environment, and aerospace [11][12][13][14][15][16]. In [17], a robust MPC strategy is presented to handle the trajectory tracking problem for an underactuated twowheeled inverted pendulum vehicle. Moreover, based on an MPC scheme, in [18], a control strategy is designed to an unmanned aerial vehicle for its automatic carrier landing system. Nevertheless, the computation complexity makes the multivariable MPC ineffectual for high speed applications where the controller must execute in a few milliseconds [19][20][21][22][23]. Moreover, the problem becomes much more complicated solving such an online constrained optimization problem by computing a numerical solver [24][25][26]. Several MPC techniques are used to overcome these problems. For instance, in [27], an explicit model predictive control moves major part of computation offline, which makes it enable to be implemented in real time for wide range of fast systems. Also, in [28], to reduce the online computational time, all the state trajectories are included in the optimal control problem as the constraints in the prediction horizon, then only a quadratic programming problem is solved. In [29], based on a mixed integer quadratic programming problem, the control input is calculated at each discrete time. In contrast to common MPC approaches, where an optimization toolbox is required, this work presents a relaxed performance index, in which the weighting matrices are computed online using the concept of Taylor series expansion and standard inverse distance weighting (IDW) functions. en, tracking performance under input-output constraints is well obtained, lighter computation load is achieved, and execution time to solve a Quadratic Program (QP) is reduced.
us, a computationally efficient constrained MPC for discrete-time statespace multivariable systems is obtained. e paper is organized as follows. Section 2 gives the preliminaries of the proposed MPC strategy. Section 3 describes the proposed relaxed cost function. Section 4 presents a tape transport system and a jet transport aircraft as study cases. Simulation results show the performance of the proposed MPC strategy and the execution time improvement compared with two well-known MPC strategies. Finally, Section 5 discusses the conclusions. Acknowledgments and the list of references finish the paper.

Model Predictive Control Based on Discrete-Time State Space Model
is section presents a brief review of MPC based on discrete-time state-space model. e original controller is proposed by Alamir in [7]. In this previous work, considering the predictions of the states, the control action is obtained through the solution of a constrained optimization problem by using a cost function with constant weighting matrices. At each sampling time, an optimal control problem is solved whose results are computationally expensive. e system dynamics is denoted by the Linear Time Invariant (LTI) State-Space Model taking the following structure: where x ∈ R n is the state vector, u ∈ R n u is the controlled input vector, y r ∈ R n r is the output vector, A ∈ R n×n is the state matrix, B ∈ R n×n u stands for the input matrix, C r ∈ R n r ×n is the output matrix, and k ∈ N denotes the sampling instant number. Hence, as in [7], from (1), the state predictions for the consecutive sampling instants are where Φ i ∈ R n×n and Ψ i ∈ R n×(Nn u ) are used to represent the N-step-ahead prediction map for Linear Time Invariant (LTI) systems in a compactness form (2). us, system (2) can be reformulated using the following vector-matrix notation: where and x(k) ∈ R N·n is the whole state trajectory of . . x n T , u(k) ∈ R N·n u concatenates the computed sequence of u � u 1 u 2 . . . u n u T , and y r (k) ∈ R N·n r is the output trajectory of y r � y 1 y 2 . . . y n r T : Now, any candidate sequence of actions u(k) has the corresponding future behavior of the system contained in the state trajectory x(k) and consequently the output trajectory y r (k).
Defining the projection matrix ⨿ (n,N) where n is the vector length, N corresponds to the prediction horizon, and i stands for the desired term, then the state, control, and output vectors at a specific instant k + i can be obtained as follows:

Development.
To find the best sequence of control action, in the present work, the value of the cost function J(u | x(k), w r (k), u p ) is defined over the prediction horizon [k, k + N] as follows: where considering a past control action value u p ∈ R n u , the last two terms are added to penalize the rate excursion of the control vector by the weighting matrix R ∈ R n u ×n u . Moreover, in the first term, the pondering matrix Q ∈ R n r ×n r penalizes the error between the output trajectory y r (k) ∈ R N·n r and the reference w r (k) ∈ R N·n r of the vector w r � w 1 w 2 . . . w n r T , which is expressed as follows: Using notation (7) and (8), (10) Moreover, substituting (1) and (2) in the first element of (10), the cost function is expressed in state-space representation: (11) us, to obtain the control law, the performance index (11) is rewritten in the following compact form: where en, the control law is obtained from (12) as follows: and applying the first action of the best sequence control action u opt (x(k)), the MPC state feedback is e MPC optimization problem (16) is solved at each sampling instant in which the measurements of the outputs and state variables are updated continuously:

Definition of the Constraints.
e systematic handling of constraints provided by predictive control strategies allows for significant improvements in performance over conventional control methodologies [30]. As in [7,31], the sequence of the future actions u cannot generally be freely chosen in R N·n u . At least, saturation constraints on the actuators have to be taken into account giving the following optimization problem P(x(k)): e cost function J(u | x(k)) is minimized considering the constraints arranged in g(u | x(k)). Nevertheless, solving general optimization problems involving a high number of decision variables and high number of constraints is generally a very difficult task. For this reason, in the present work, based on the Inverse Distance Weighting (IDW) method [24,[32][33][34][35][36], the cost function (12) input and output constraints are included in Q and R matrices, where, to reduce computational complexity, Taylor series expansion is used to obtain the predicted manipulation u p and predicted output y p [23,[37][38][39][40][41]: e penalization at the output Q ∈ R n r ×n r is defined: where q p � e σ · c, ∀p ∈ 1, . . . , n r , α is the tuning parameter set for a desired closed-loop performance, σ � α + |y p − w p |min(|y p − y max p |, |y p − y min p |), and c � 1/(min(|u p − u max p |, |u p − u min p |)), which corresponds to an IDW function.
Mathematical Problems in Engineering e deviation between the predicted output y p and the reference w p is penalized by a function of distance |y p − w p |. Furthermore, to avoid that the predicted process variable y p remains close to the lower bound y min p or the upper bound y max p the min(|y p − y max p |, |y p − y min p |) term is added. Finally, c term is included in order to increase the penalization while the predicted manipulation variable u p is close to the upper or lower bounds, u max p and u min p , respectively. e penalization at the input R ∈ R n u ×n u is defined: where r q � c + |δu q − δu p q |, ∀q ∈ 1, . . . , n u . e first term is included to avoid that the predicted manipulated variable u p remains close to the upper or lower bounds, u max q and u min q , while the second term is added to reduce the large excursion between the increment of the predicted manipulated variable δu q and the increment of the previous manipulated variable δu p q . e resulting optimization problem P(x(k)) expressed as the cost function J(u | x(k), w r (k), u p ) is minimized at each instant k by a sequence of future actions u opt (x(k)) respecting all the constraints. As it is described, Q and R are defined as functionals of constraints depending on the current state x(k), the desired trajectory w r (k), and the past controlled input u p . us, the matrices H, F i , i � 1, 2, 3 are computed online looking for relaxing the optimization problem.

Algorithm for the Proposed MPC Strategy.
e algorithm used for the proposed discrete-time state-space constrained MPC strategy consists of (1) Define the linear time invariant mathematical model of the physical system in state-space representation (equation (1)) (2) Compute the matrices Φ i and Ψ i , ∀i ∈ 1, . . . , N { } used to represent the N-step-ahead prediction map for LTI systems (equation (4)) (3) Estimate the manipulation u p (equation (18)) and the output y p (equation (19)) at the next sampling instant, using Taylor series expansion (4) Compute the weighting matrices Q (equation (20)) and R (equation (19)) based on the proposed functions (5) Compute the online matrices H, F i , i � 1, 2, 3 of the proposed performance index (equation (12)) (6) Obtain the control law u opt (x(k)) (equation (14)) and apply to the system the first action of the best sequence control action K MPC (x(k)) (equation (15)) (7) At the next sampling instant (k ≔ k + 1 ), the reference w r , the bounds y min , y max , u min , and u max , the measurements of the outputs y r , and the state variables x(k) are updated, and the MPC optimization problem is solved again (equation (17)); the algorithm goes to step 3

Simulation and Results
e present work and the MPC strategies in [7,31] are developed using the same computational platform for its evaluation. Hence, to solve the optimization problem in [7,31], quadprog toolbox from MATLAB ® is used. MATLAB ® codes for the following study cases are available: MPC Matlab Files_MPE.

Example 1: Tape Transport
System. e tape drive system consists of two reels to supply and file data. Here, the data transfer rate is proportional to the tape transport speed.
us, the tape drive mechanism must be able to rapidly transport a fragile tape with an accurate tension regulation. Figure 1 shows the schematic of a tape transport system where its components and variables involved are the tape stiffness and the damping denoted by K and D, the reel radii and the inertia represented as r and J, the motor torque constant K t , and the viscous friction coefficient denoted by β.
Assuming there is no force loss across the head, the tape tension T � T 1 � T 2 [42]. Although the physical model of the process contains nonlinearities, in (22), a simplified state-space model is presented in continuous time [42][43][44]: e model has three states, x � T ω 1 ω 2 T , where T is the tension tape in N; meanwhile, ω 1 and ω 2 represent the 4 Mathematical Problems in Engineering supply and take-up reel in rad/s, respectively. Moreover, the system has two inputs u � u 1 u 2 T that represent the voltages applied to the reel motors in volts, and two outputs y � y 1 y 2 T which stand for the tape speed v rw at the readwrite head in m/s and the tape tension T, respectively. e control strategy described in the present work is simulated using parameters from the tested tape system described in [42,45,46], whose parameters are summarized in Table 1.
Considering that the motors are nominally identical, for both motors, it is used as the same motor torque constant K t and viscous friction coefficient β [45,46]. Discretizing (23) with a sampling time τ � 0.1 s and considering a prediction horizon N � 10 with a tuning parameter α � 3.7, the simulation results are shown in Figure 2. It is divided in two main parts, the first part corresponds to the outputs and the inputs of the system using the MPC with the relaxed cost function, while the second part shows the coefficient values of the weighting matrices Q and R. As it is shown, computing the coefficient values q 1 and q 2 , the constrained tape speed y 1 and tape tension y 2 present a good performance. Here, y 1 do not present overshoot and has a maximum settling time of 0.1 seconds while y 2 has a maximum overshoot of 3% and a settling time of 0.3 seconds. On the contrary, the motor voltages u 1 for the supply reel and u 2 for the take-up reel remains inside their lower and upper bounds by using the computed coefficient values r 1 and r 2 . Here, the values of r 1 are lower than the values of r 2 , due to the rates of the manipulations.
In order to see the closed-loop stability of the system, the stability indicator S N is defined as in [7]: where λ j stands for the eigenvalues of the system for j ∈ 1, . . . , n is the closed-loop gain. Figure 3 shows that the system remains stable during the test. Finally, Table 2 is presented to compare the execution time between the present work and previous works [7,31].
As it can be seen, the total execution time is reduced, by taking advantage of the relaxed cost function. Here, the computation of u(k) is 1.924 seconds and 2.148 seconds faster than the computation of the manipulations using [7,31] MPCs strategies. Henceforth, the percentage consumption of time to obtain the control actions is notably decreased.

Example 2: Jet Transport Aircraft.
e Jet Transport Aircraft Boeing 747 in high-lift configuration addresses complex geometries and physical phenomena that make the controller design a difficult process. Figure 4 illustrates the Jet Transport Aircraft with its components and variables involved such as the angles β and φ and the angular velocities ψ and θ.
Although the physical model of the Boeing 747 is lengthy, in (24), the simplified state-space model during cruise flight at MATCH � 0.8 and H � 40, 000 ft is presented in continuous time [47]: e model has four states, x � β ψ θ φ T , where β is the sideslip angle, φ stands for the bank angle, and meanwhile ψ and θ represent the yaw and roll rate, respectively. Herein, all the angles are in rad and the angular velocities in rad/s. e system has two inputs u � u 1 u 2 T : the rudder and the aileron deflections, and two outputs y � y 1 y 2 T : the yaw rate ψ and the bank angle φ.
Using a sampling time τ � 0.2 s, system (24) is discretized. en, a set of changes in the output reference and a series of variations in the constraints are used to test the present work. Considering a prediction horizon N � 20 and a tuning parameter α � 4.5, simulation results are shown in Figure 5. Here, using the computed coefficient values q 1 and Supply reel Take-up reel  q 2 , the yaw rate y 1 has a maximum overshoot of 2.9% with a maximum settling time of 1 second, while the bank angle y 2 has a maximum overshoot of 6.6% with a maximum settling time of 1.2 seconds. Moreover, with the computed coefficient values r 1 and r 2 , the rudder deflection u 1 , and the aileron deflection u 2 remains inside their bounds without saturation. Here, the values of r 2 are greater than the values of r 1 , this is because the rate of the manipulation u 2 is greater than the rate of the manipulation u 1 . On the contrary, Figure 6 shows the stability behavior during the test. As it is shown, the system remains stable. Finally, the execution time comparison between the present work and previous works [7,31] is shown in Table 3. As it can be seen, the total execution time is considerably reduced due to the relaxed cost function. Here, the computation of u(k) is 6.074 seconds and 3.155         Mathematical Problems in Engineering seconds faster than the computation of the manipulations using [7,31] MPCs strategies. en, as in example 1, the percentage consumption of time to obtain the control actions is notably decreased.

Conclusions
is paper presents a discrete-time state-space MPC approach for multivariable systems. Based on the IDW method and the concept of Taylor series expansion, a relaxed performance index with constraints defined in the online weighting matrices is proposed to compute the control action. us, as in study cases, the proposed MPC strategy is used to control a tape transport system and a jet transport aircraft during cruise flight.
Simulation results show that the proposed MPC strategy with the relaxed cost function has a good performance, no matter abrupt changes of set-points and constraints occur, even at the same time. Additionally, compared with two well-known discrete-time state-space MPC strategies, there is a significant improvement on the execution time without affecting the tracking performance. e percentage consumption of time to compute the best sequence of control actions u(k) is 1.3% for the tape transport system and 1.7% for the jet transport aircraft. Henceforth, it takes almost 0.1 milliseconds for the tape transport system and 0.12 milliseconds for the jet transport aircraft to obtain the manipulation u(k) that minimizes the proposed cost function while respecting the constraints. us, the proposed MPC strategy with the relaxed cost function is mainly suitable for constrained multivariable real processes with fast dynamics.
Data Availability e data of the conducted experiments and simulations are available upon requirement.

Conflicts of Interest
e authors declare that they have no conflicts of interest.