Abstract

Model Predictive Control (MPC) can effectively handle control problem with disturbances, multicontrol variables, and complex constraints and is widely used in various control systems. In MPC, the control input at each time step is obtained by solving an online optimization problem, which will cause a time delay in real time on embedded computers with limited computational resources. In this paper, we utilize adaptive Alternating Direction Method of Multipliers (a-ADMM) to accelerate the solution of MPC. This method adaptively adjusts penalty parameter to balance the value of primal residual and dual residual. The performance of this approach is profiled via the control of a quadcopter with 12 states and 4 controls and prediction horizon ranging from 10 to 40. The simulation results demonstrate that the MPC based on a-ADMM has a significant improvement in real-time and convergence performance and thus is more suitable for solving large-scale optimal control problems.

1. Introduction

Model Predictive Control (MPC) [1, 2], a.k.a. receding horizon control, is a classical optimization-based control strategy for multivariable constrained systems. Since its appearance in the 1970s, MPC has been widely applied in a variety of areas, such as the transportation system [3, 4], grid system [5], and water distribution system [6]. MPC only utilizes the first element of the control input sequence obtained by solving an optimization problem. On the one hand, although receding horizon way makes MPC robust to model bias and external noises, MPC requires solving the optimization problem in each control cycle, which may cause a large time delay. On the other hand, low-complexity optimization algorithms have been motivated by the fact that several of these recent applications exploit microcontroller that has limited computing power [79].

For linear systems with a quadratic cost function and multiple linear constraints, MPC can be easily transformed into a standard quadratic programming (QP) form [10]. Therefore, off-the-shelf QP solvers can be used in this case. Most of them can be categorized into one of the two main classes, i.e., active-set and interior-point methods. Interior-point methods have become common choice for implementation because of its good convergence performance [11, 12]. However, it tends to demand significant online computation effort to compute the input. Recently, some modified methods have been proposed to reduce the computational complexity [13]. For example, Domahidi et al. decreased the computational burden by introducing a low-rank matrix forward substitution scheme and extended the interior-point methods to applications on low-cost embedded hardware. Compared to interior-point methods, active-set methods need on average substantially more iterations, but each iteration is computationally much cheaper and can be ward- or hot-start. For example, a parametric active-set algorithm for quadratic programming (qpOASES) is another commonly used QP solver which is applicable when numerous related QPs are solved sequentially, by exploiting the geometric property of state space [14]. Unfortunately, since qpOASES is an active-set-based method, it is limited to small- to medium-scale QP problem and sensitive to the external disturbance [14].

Aforementioned methods are centralized optimization algorithms, which are not appropriate for real-time large-scale applications on low-cost embedded computers. On the contrary, due to its success in solving large-scale optimization problems and pretty low computational complexity, distributed optimization algorithms have drawn increasing attention [15]. As a distributed optimization algorithm, Alternating Direction Method of Multipliers (ADMM) is a strong algorithm for solving convex optimization problems [16, 17]. ADMM was first introduced in the mid-1970s for the numerical solution of partial differential equations [18, 19] and particularly useful for solving optimization problems that are too large to be handled by generic optimization solvers. Because of its fast processing speed and good convergence performance, the method has been used in a large number of areas, such as compressed sensing [20], neural network [21], and image processing [22]. More details can be found in [17]. This wide range of applications has triggered a great interest in developing a better and more comprehensive understanding of the theoretical properties of ADMM.

ADMM has been proved to be an excellent fit for QP problems, since it can achieve linear convergence with strong convexity assumptions [23]. However, a drawback of ADMM is that the number of iterations required to converge is highly dependent on the step-size parameters. Unfortunately, analytic results for the optimal choice of penalty parameter are not available except for very special problems [24, 25]. In the special case of a strictly convex quadratic problem, criteria for choosing an optimal constant penalty have been recently proposed in [24, 26]. Although the optimal penalty parameter can significantly improve practical convergence, these methods make strong assumptions about the objective and constraint matrix. By introducing relaxation parameter, a method called relaxed ADMM uses past iterations to compute the current one. Instead of the fixed penalty parameter, Aybat and Iyengar [27] suggested using a penalty parameter sequence. Another different approach to the above methods is to accelerate ADMM by adaptively modifying the penalty parameter. As one of the adaptive ADMM methods, residual balancing (RB) [28] tunes penalty parameter to keep both residuals of similar magnitude. We borrow the idea from this method to further improve the convergence of ADMM combined with MPC.

In this paper, we utilize a-ADMM to efficiently solve the MPC problem, and the performance is verified by multiple simulation experiments. To be more specific, we transform the standard MPC problem into a distributed formulation and then take advantage of a-ADMM to solve this problem. We verified the method on an Unmanned Aerial Vehicle (UAV) control problem, and this mechanical control system includes 12 states and 4 control inputs [29]. We tested our method on different settings, i.e., the prediction horizon N ranges from 10 to 40. For each setting, MPC requires the solution of a QP with variables and inequality constraints. In this complex UAV control problem, a-ADMM algorithm achieves an absolute improvement of real-time response ability over conventional interior-point method and convergence performance over other ADMM variations.

The paper is organized as follows. Section 2 briefly introduces MPC, ADMM, and some ADMM variations. Section 3 describes the main procedures of the transformation from MPC to the distributed optimization problem and the application of a-ADMM in this specific problem. In Section 4, we compared the method with some baselines in the UAV control task. Finally, the conclusion is drawn and our future work is declared in Section 5. Appendix A describes how to formulate the QP problem of MPC in some details.

2. Backgrounds

2.1. Model Predictive Control

MPC is an optimal control approach which can handle multiple constraints on states and control inputs. Using a receding horizon control framework, MPC is robust to the model bias and external noises, which makes MPC fit well for highly complex multivariable industry processes.

In this paper, the following discrete-time linear time-invariant system is considered:where and denote the state and control input, respectively, at time step t. and are the system and control matrices, respectively. For the regulation problem, at each time step t, MPC solves the finite horizon optimal control problem:where N is the prediction horizon. denotes the predicted state vector at time , obtained by applying the input sequence . and define constraints on state variable , and define constraints on input variable . and are stage and terminal matrices which penalize the state deviation at the end of the prediction horizon N and over the entire horizon, respectively. is the cost matrix for the control inputs. Typically, , , and .

In MPC, control inputs are obtained by solving open-loop optimal control problem (2a)–(2e) [10]. To be more specific, at each time step t, the optimization problem (2a)–(2e) is solved by taking the current state as initial state of the problem and generates a sequence of optimal control input. However, according to receding horizon framework of MPC, only the first element of the control input sequence is applied to the plant. The remaining solutions are culled, and optimal control problem with new initial state will be solved again at the next time step .

2.2. Alternating Direction Method of Multipliers

ADMM is a strong algorithm which can efficiently solve convex optimization problems, and it can decompose the global problem into several smaller and easy-to-solve local subproblems [30]. Because of its fast processing speed and good convergence performance, ADMM has received tremendous interest and approval for solving numerous problems in machine learning, statistics, and signal processing [17, 3133].

ADMM can efficiently solve the optimization problems of the following standard form:where f and C are convex. This problem can be rewritten in another form aswhere is the indicator function of C, is primal variable, and is separate variable. The augmented Lagrangian is defined aswhere ρ is the step size (or penalty parameter) and c is the scaled dual variable. The general augmented Lagrange multiplier method is to minimize , w.r.t. and at the same time. This is in difficulty and does not exploit the fact that the objective function is separable. To remedy this issue, ADMM decomposes the global problem (4) into two subproblems and minimizes and , respectively, with the following three steps:

Under rather mild conditions, ADMM converges for any constant penalty parameter [34]. The convergence of ADMM is based on the primal and dual residuals:

It has been observed that the algorithm approaches to the optimal solution when the primal and dual residuals approach to zero. Generally, the iterative process is terminated ifwhere is absolute tolerance and is relative tolerance. However, ρ has an immediate effect on the convergence element of the algorithm, and inappropriate tuning of this algorithm parameter may render the method moderative.

One way of accelerating the convergence properties of the algorithm is to utilize previous iterates when computing the following ones. This practical method is called relaxed Alternating Direction Method of Multipliers (r-ADMM) [35] that replaces with in the and updates (6a)–(6c), yielding

The parameter is called the relaxation parameter. Note that letting for all k recovers the original ADMM iterations. Empirical studies have illustrated that overrelaxation, i.e., letting , is often beneficial and the guideline has been proposed [36].

Another approach for improving the convergence properties of the algorithm is to utilize different penalty parameters for each iteration. Adaptive Alternating Direction Method of Multipliers (a-ADMM) [28] is based on the following observation: increasing strengthens the penalty term and leads to smaller primal residuals but larger dual ones; on the contrary, decreasing results in smaller dual residuals and larger primal residuals. As both residuals must be small at convergence, it requires adjusting ρ to keep both primal residuals and dual residuals of the same magnitude. A simple and effective scheme for this purpose is as follows:where μ > 1, >1, and >1. Using different penalty parameter for each iteration k, a-ADMM can converge quickly and be insensitive to the initial penalty parameter. However, without a careful choice of μ, , and , this algorithm may fail to converge [28].

3. Fast MPC Based on Adaptive ADMM

3.1. Adaptive ADMM MPC

We focus on improving MPC by solving its QP problem with a-ADMM. In this work, the standard MPC problem is transformed into an equivalent QP problem which could be easily decomposed into several smaller local subproblems. And a-ADMM could be finally adopted to efficiently solve the resulting problem.

MPC problem ((2a)–(2e)) can be expressed as a QP problem when the system dynamic and constraints are linear and the cost function is quadratic [10]. Specifically, by substituting into the cost function, the optimization problem equations (2a)–(2e) can be formulated aswhere is the optimization vector and H, C, Y, , , and S are constant matrices that can be easily obtained offline. Note that H has to be positive definite for the convexity of this problem. The detailed process of this transformation is shown in Appendix A.

As the problem is determined by the current state , the implementation of MPC requires the online solution of problem ((11a) and (11b)) at each time step. Although off-the-shelf QP solvers, e.g., active-set methods [37] and interior-point methods [11], are available, it may demand significant computation effort in order to compute the input online.

Fortunately, the QP problem (11a) and (11b) can be further transformed into ADMM standard form by introducing a split vector named slack vector (see [17]):where is the initial objective function with qualified domain and is defined as the indicator function:

The augmented Lagrangian is defined aswhere ρ is the penalty parameter and is the dual parameter. We propose to solve problem (11a) and (11b) with a-ADMM as follows:where the penalty parameter ρ is a piecewise linear function as equation (10). As the optimizer is only needed, the term involving Y is usually removed from problem (11a) and (11b). Since both residuals must be small at convergence, ρ is tuned to keep both residuals of similar magnitude.

3.2. Convergence Analysis of Adaptive ADMM-MPC

He et al. [28] proved that convergence is guaranteed for a-ADMM when either of the two following conditions is satisfied:Condition 1 (bounded increasing):where .Condition 2 (bounded decreasing):where .

Condition 1 (Condition 2) suggests that the increasement (decreasement) of adaptive penalty parameter is bounded.

Based on the above analysis, the framework of a-ADMM algorithm is portrayed in Algorithm 1.

a-ADMM algorithm
Initialization: optimization variables , ,
While not converge by stopping criteria ((8a) and (8b)) and k¡maxiter do
if
  
else if
  
else
end if
end while
Output: optimal control input

4. Simulation Experiment

We evaluated our approach on an UAV control task. The dynamic model of UAV used in this experiment is from [29] and can be expressed as follows:where x, y, and z denote the position of the UAV and θ, ϕ, and ψ denote the rotation angle of the UAV around the Cartesian coordinate axis. The damping coefficient β takes into account the actual friction effect. This dynamic model has 12 states and 4 control inputs . For the purpose of designing a appropriate MPC controller for UAV, the nonlinear dynamic model ((18a)–(18f)) is linearized at the equilibrium point (see [29]):where is the state vector, is the input vector, and is the output vector. Assuming that the UAV system is completely measurable, and the observation matrix C is a 12 × 12 identity matrix. The coefficient matrices A and B are as follows:

The baseline methods we compared in our experiments are as follows.

4.1. Interior-Point Method (IPM)

IPM is the mostly used QP solver modeling the problem constraints as parametrized penalty functions which also referred to as barrier functions. At each iteration, an unconstrained optimization problem is solved for varying barrier function until the optimum is achieved.

4.2. qpOASES: a Parametric Active-Set Algorithm for Quadratic Programming

qpOASES is a very used QP solver which is applicable when numerous related QPs are solved sequentially, by exploiting the geometric property of state space [14].

4.3. Alternating Direction Method of Multipliers (ADMM)

ADMM is a strong algorithm for solving convex optimization problems [16, 17] which decomposes the global problem into several smaller and easy-to-solve local subproblems [30].

4.4. Relaxed Alternating Direction Method of Multipliers (r-ADMM)

r-ADMM is one of the variants of ADMM which introduces relaxation parameter α to relax the update of the optimal solution. Empirical studies have illustrated that overrelaxation is beneficial to the convergence [36].

We first fixed the prediction horizon to find the optimal setting of ADMM and r-ADMM: the penalty parameter ρ and the relaxation parameter α. Figure 1 gives the results of different values of parameters. As can be seen, the average number of iterations depends significantly on the penalty parameter ρ and the relaxation parameter α. Based on this result, we simply set and in the following experiments.

Similarly, we used the same control task to find the optimal setting of a-ADMM. According to the theory of a-ADMM, the relaxation parameter α is fixed and can be set α = 1.5 based on previous results. However, different from ADMM, the penalty parameter ρ of a-ADMM is a piecewise linear function as equation (10) influenced by three other parameters, i.e., μ, , and . Figure 2 gives the results of different values of μ, , and on the UAV control task when the prediction horizon . Based on the lowest average iterations of each setting, we set , , and in the following experiments, and the optimization variables , , and are all initialized to a zero vector, so that the initial strategy of optimization variables does not affect the performance comparison of these algorithms (IPM, qpOASES, ADMM, r-ADMM, and a-ADMM). The parameters of these algorithms are listed in Table 1.

The feedback loop controls the quadcopter to land on a platform located at coordinate position . As shown in Figure 3, the closed-loop performance (including the feasibility and stability) of ADMM and a-ADMM is exactly the same. The reason is that, with the same absolute tolerance and relative tolerance , ADMM and a-ADMM converge to the same optimal solution leading to the same closed-loop performance.

Figure 4 gives our real-time performance of each method on the UAV control task. This figure shows the cumulative running time of our method and other baselines with the prediction horizon N ranging from 10 to 40. From the simulation results, it can be seen that, combined with MPC, ADMM and its variants significantly outperform custom QP solver, IPM, and qpOASES. With the increment of the task complexity, the running time of IPM and qpOASES increases obviously, while that of ADMM, r-ADMM, and a-ADMM basically remains unchanged.

As shown in Figure 4, there is not much difference in the cumulative running time of several variants of the ADMM algorithm. In this part, we aim to further investigate the difference in the cumulative number of iterations between them. Figure 5 portrays the cumulative iteration number of each method when the prediction horizon N ranging from 10 to 40. Apparently, although ADMM and its variations have the similar real-time performance, a-ADMM significantly outperforms all other baselines in terms of convergence performance under all the task settings. qpOASES proposes to move on a straight line in the state space when transitioning from the “old” QP problem to the “new” one. The reason why the cumulative iteration number of qpOASES stays stable after several time steps is mainly because the state of quadcopter has converged to the objective region, and the adjacent QP problem is quite similar.

5. Conclusions

In this paper, we utilize a-ADMM to accelerate the solution of MPC. We firstly transform the standard MPC problem into a distributed formulation and then take advantage of a-ADMM to solve this problem. We verified this method on an UAV control problem and compared a-ADMM with IPM, ADMM, and r-ADMM under different control settings, i.e., the prediction horizon N ranging from 10 to 40. The experimental results have shown that a-ADMM can greatly improve the effectiveness of MPC. As a future and noteworthy research direction, one might provide better scheme for adaptively adjusting penalty parameter considering more complex constraints.

Appendix

A. Equivalent Representation of

The linear MPC problem can be expressed as QP problem by substituting into cost function. Firstly, we verified that has an equivalent representation based on the following transformation:

The conclusion can be deduced, . Having disposed of this preliminary step, we begin to prove the equivalent form of the objective function in Appendix B.

B. The Transformation of the Objective Function

Secondly, we transform the objective function into its equivalent form. As previously mentioned, the objective function can be casted into by substituting , and the formula is derived as follows:

Set and , and can be expressed as consequently,where

In a similar way, can be rewritten and finally obtained aswhere

C. The Transformation of Constraints

The last part is about the transformation of constraints, and the formula is derived as follows:

By substituting , the constraints are equivalent to through the following process:

By substituting formula (A.8) into formula (A.9), we obtain

The solution of is obtained as follows:

That is,

Formula (A.12) is equivalent towhere

Data Availability

We use simulation data, and our model and related hyperparameter are provided in our paper.

Conflicts of Interest

The authors declare that they have no conflicts of interest.

Acknowledgments

This work was jointly supported by the National Natural Science Foundation of China (11871077), the Natural Science Foundation of Anhui Province (1808085MA04), and the Natural Science Foundation of Department of Education of Anhui Province (KJ2017A362).