An integrated path-tracking and control allocation method for autonomous racing electric vehicles

Abstract

In recent years, path-tracking controllers for autonomous passenger vehicles and Control Allocation (CA) methods for handling and stability control have both received extensive discussion in the literature. However, the integration of the path-tracking control with CA methods for autonomous racing vehicles has not attracted much attention. In this study, we design an integrated path-tracking and CA method for a prototype autonomous racing electric vehicle with a particular focus on the maximising the turning speed in tight cornering. The proposed control strategy has a hierarchical structure to improve the computational efficiency: the high-level path-tracking Model Predictive Control (MPC) based on a rigid body model is designed to determine the virtual control forces according to the desired path and desired maximum velocity profile, while the low-level CA method uses a Quadratically Constrained Quadratic Programming (QCQP) formulation to distribute the individual control actuator according to the desired virtual control values. The proposed controller is validated in a high-fidelity simulation vehicle model with the computational time of the optimisation controller presented to demonstrate the real-time control performance.

KEYWORDS:

• Autonomous vehicle
• control allocation
• MPC
• path-tracking

1. Introduction

The autonomous vehicle has been considered as the next generation vehicle and has received increased attention in the recent literature. A vital part of the control architecture of an autonomous vehicle is around the path following control, which includes two parts: the high-level motion planner, which uses environmental sensing and perception techniques to determine a desired path and speed profile, and the low-level path-tracking controller which is designed to track the desired path with desired speed and maintain the vehicle dynamics stability. In this study, we focus on the low-level path-tracking controller.

Many of earlier studies about the path-tracking used a driver model to determine the steering control input [1–5], then the study on path-tracking for driver model has been switched to autonomous passenger vehicle application recently [6–11]. In order to track the desired path and minimise the lateral tracking error, Attia et al. utilised a Nonlinear Model Predictive Control (NMPC) to determine the control input of steering angle [6]. In another study, based on the lateral tracking error, the H-infinity technique was applied to design the lateral tracking controller [7]. Tagne et al. proposed a higher-order Sliding Mode Controller (SMC) for the reference path-tracking control, since the robustness of SMC can be beneficial to overcome the nonlinear terms of the vehicle model [9].

In recent years, the independent wheel control in an over-actuated vehicle has been widely discussed in the vehicle dynamics and control literature. The independent wheel control provides redundant control actuators and is becoming popular due to the possibility of enhancing the performance of existing controller through the carefully design of Control Allocation (CA) strategy [12]. It is also suggested that CA can improve the cornering response by shaping the under-steer characteristics and ensuring stability in limited-handing condition [13]. Thus, integrating the CA strategy with path-tracking controller can greatly improve path-tracking and vehicle dynamics performance of autonomous vehicle with in-wheel motor, which has become a very attracting research area.

CA was originally used to control and allocate traction or brake toque of individual tire for the vehicle handling and stability control of traditional human-driving passenger vehicles [14–18]. In recent years, quite a few approaches integrate path-tracking with CA based on Direct Yaw Control (DYC) for autonomous passenger vehicles [19,20]. Guo et al. proposed a Linear Time Varying (LTV) based MPC as the path-tracking controller to determine the steering angle and required yaw moment, then a low-level CA law is applied to distribute the individual wheel torque to achieve the required yaw moment [19]. Wang et al. developed a modified Composite Nonlinear Feedback (CNF) controller to track the desired yaw rate and lateral velocity and then the required yaw moment control value is distributed to longitudinal tire force of individual wheel [20]. In [13], a hierarchical structure and a holistic structure of path-tracking controller is suggested and integrated with torque vectoring controller. In the above studies, a distribution of the DYC command on the individual wheel torque (also called torque vectoring) is applied. However, the CA problem can be extended to the distribution of the steering angle, motor driving torque, motor regenerating brake torque and hydraulic brake torque on the individual wheel. In this way, the redundant vehicle control actuators can be fully utilised to achieve the best control performance, but CA problem becomes more complex and requires a carefully designed CA strategy. In our previous study [21], the hierarchical control structure is proposed by including high-level PID based path-tracking controller and low-level CA approach. The CA method can optimally distribute the individual wheel torque and front wheel steering angle. Nevertheless, the simple PID path-tracking controller can only track the desired path with fixed longitudinal speed, which cannot optimise the path-tracking with maximum cornering speed in racing scenario (without the help of MPC). Furthermore, the CA method still has not fully considered the control allocation of front and rear wheel steering angle, motor driving torque, motor regenerating brake torque and hydraulic brake torque of individual wheel.

MPC has been widely applied in the path-tracking control of autonomous passenger vehicles since the control performance can be greatly improved by MPC through the preview of the future road profile [19,22–27]. However, there are limited studies on the path-tracking control of the autonomous racing vehicles. Alcala et al. proposed a linear parameter varying (LPV) MPC for autonomous racing vehicles [28]. The validation of the controller is based on the scaled electric vehicle and the vehicle speed is relatively low. It is questionable if the proposed controller can achieve good control performance for the high-speed full-sized autonomous racing vehicle. Similarly, Liniger et al. proposed a linearised MPC for the path-tracking of autonomous racing vehicle and an integrated MPC contouring control to achieve the combined path planning and path-tracking control, where the nonlinear bicycle model considering the longitudinal and lateral friction forces is modified as the LTV model [29]. In order to improve the internal model accuracy of the MPC, the learning based MPC approaches are proposed in [30]. The internal MPC model includes a nominal simple vehicle model and a self-learning model to improve the vehicle model accuracy. The major issue of this data-driven self-learning approach is the uncertain stochastic internal model will render the predicted state a random variable. This uncertainty will propagate throughout the whole prediction horizon, and the prediction model system may be unstable with a longer prediction horizon. Furthermore, the Gaussian Process Regression has been applied to train the self-learning model and the large number of training data points will significantly increase the computational burden.

In summary, in the current literature, the path-tracking controller design for autonomous racing vehicles is mainly based on the linearised MPC approach with simple internal model to improve the computational efficiency and guarantee the controller reliability. Although the self-learning MPC approach has been tried to improve the internal model accuracy, the reliability of this approach is questionable. Although CA and MPC approaches have been extensively applied in the path-tracking control for passenger vehicles, the integration of the path-tracking control and maximum speed control in tight cornering for the racing vehicle is less focused. The optimisation of redundant steering and driving actuators can greatly enlarge the working envelope of racing vehicles and push the racing vehicles close to the friction limit.

The aim of our study is to design an integrated path-tracking and CA method for a prototype autonomous racing electric vehicle in limited handling scenario, with the proposed method being able to satisfy the requirement for real-time implementation. The prototype vehicle under consideration has independent front- and rear-wheel steering, two hub motors with driving/regenerating braking capabilities on the two rear wheels, and one motor on the front axle that through an open differential equally distributes the driving/regenerating torque among the two front wheels. In addition, there are two independent hydraulic brake channels, one on the front axle and one on the rear axle of the vehicle, which is shown in Figure 1. Motivated by the above shortcomings in current literature, a hierarchical control structure of integrated path-tracking and CA method is presented in this study. In the high-level controller, a path-tracking MPC is proposed based on a rigid body vehicle model including the virtual control inputs of the total longitudinal, lateral force and yaw moment on the vehicle, while individual tire friction force and tire model are considered in the low-level controller and decoupled from high-level to improve the computational efficiency. This way a modular strategy is achieved, which can be easily modified and applied to most of the existing actuator structures of an electric vehicle.

Figure 1. The structure diagram of the complete control system.

The major contribution of this study can be summarised as follows. First, the proposed hierarchical controller applies the virtual forces and moment to performance high-level path-tracking motion and effectively move the effect of the actual friction tire forces in the lower CA level which also accounts for the various actuators. Secondly, the path-tracking MPC with prediction horizon is proposed based on a rigid body vehicle dynamics model, and the tire friction circle and maximum yaw moment constraints are included in the high-level MPC to guarantee the feasibility of the computed commands to be passed on the low-level CA. Furthermore, a two-layer CA method is proposed to control and allocate the various steering and driving actuators: the upper layer CA optimally distributes the individual longitudinal tire force and steering angle according to the desired virtual control values from the path-tracking controller, formulated as a Quadratically Constrained Quadratic Programming (QCQP) problem to improve the computational efficiency, while the lower layer rules-based method is designed to distribute the required individual longitudinal tire force into the individual motor and hydraulic brake torques. Finally, the time-efficient solver generated by Forces Pro is used to solve the high-level path-tracking MPC and optimisation problem of CA in low-level in this study, which can satisfy the requirement for real-time path-tracking control [31]. It is noted that our study only focuses on the high-level path-tracking and low-level CA algorithm, and the detailed actuator model of the powertrain system is out of the scope and not considered in our study. The major novelty of our study can be suggested: (1) the path-tracking control and maximum speed control in tight cornering (through the offline pre-calculation of maximum cornering speed and the online speed control in MPC prediction horizon) are integrated together for the autonomous racing vehicle; (2) through more powerful in-wheel motor system and the potent racing tire (tire-road friction coefficient > 1), the proposed controller can fully utilise the potential of the racing vehicle and reach the acceleration limit; (3) in the high-level MPC path-tracking controller, the virtual total desired tyre force and yaw moment (yaw moment constraint is not well addressed in the literature) are constrained so that the low-level CA can successfully allocate individual actuators.

This paper is organised as follows: the high-level path-tracking MPC is presented in section II, then low-level CA strategy is introduced in section III, and finally simulation examples are used to validate the effectiveness of the proposed integrated control method in section IV.

2. High-level path-tracking MPC design

The proposed hierarchical controller includes the high-level path-tracking controller and low-level CA method. The high-level path-tracking MPC method tracks the desired path and velocity profile and sends the required virtual control forces into the low-level CA to distribute the individual control actuator value to the vehicle. Vehicle measurements are also used to feed back the necessary vehicle state values to the high-level and low-level controllers. The detailed structure of the complete control system is presented in Figure 2.

Figure 2. The structure diagram of the complete control system.

2.1 Vehicle dynamics model for high-level path-tracking MPC

In this section, a three degrees-of-freedom (DOF) rigid-body nonlinear dynamics model considering the combined longitudinal, lateral and yaw dynamics is used to design the path-tracking MPC. The virtual control inputs of total longitudinal tire force, lateral tire force and yaw moment are utilised to simplify the vehicle dynamics model. In addition, the yaw angle error and lateral tracking error are also included in the vehicle model to present the path-tracking performance [20]. The 3 DOF vehicle model is presented by the following equations: (1a) $\begin{array}{rl}{\dot{v}}_x& =v_{y}r+\frac{F_{\mathrm{xd}}}{m}\end{array}$(1a) (1b) $\begin{array}{rl}{\dot{v}}_y& =-v_{x}r+\frac{F_{\mathrm{yd}}}{m}\end{array}$(1b) (1c) $\begin{array}{rl}\dot{r}& =\frac{M_{\mathrm{zd}}}{I_z}\end{array}$(1c) (1d) $\begin{array}{rl}{\dot{\psi }}_e& =r-\kappa v_x\end{array}$(1d) (1e) $\begin{array}{rl}{\dot{Y}}_e& =v_x\sin {\psi }_e+v_y\cos {\psi }_e\end{array}$(1e) where $v_x,v_y,r$ represent the vehicle longitudinal velocity, lateral velocity, and yaw rate, respectively, $I_z$ and $m$ are the moment of vehicle inertia in terms of yaw axis and vehicle mass, $F_{\mathrm{xd}},F_{\mathrm{yd}},M_{\mathrm{zd}}$ are total longitudinal tire force, lateral tire force and yaw moment, $Y_e$ and ${\psi }_e$ are the vehicle lateral deviation and yaw angle error from the reference path and $\kappa$ is the curvature of desired path.

2.2 The design of high-level path-tracking MPC

The general formulation of a nonlinear MPC can be presented as following equations by selecting the best control input to satisfy the cost function $J$: $\underset{x,u}{min}\sum _{k=0}^{N-1}J(x_k,u_k,p_k)$Subject to: $x_{k+1}=f(x_k,u_k,p_k)$ (2a) $\begin{array}{rl}& k=0,1,\dots N-1\end{array}$(2a) (2b) $\begin{array}{rl}& x_k\in \mathrm{Xk}=1,\dots N-1\end{array}$(2b) (2c) $\begin{array}{rl}& x_0=x(t)\end{array}$(2c) where $N$ is the horizontal number of predicted time steps, $x_k$ is the predicted vehicle state in the prediction horizon and $u_k$ is the predicted control input, $f(x_k,u_k,p)$ in equation (2a) is the nonlinear system prediction model, which is the discrete formulation of nonlinear dynamics system (1) and is formulated in (4) and $p$ is the external changing parameters determined outside of MPC. The boundary constraints of vehicle states are presented in (2b) while the initial values $x_0$ can be determined by the measured feedback state values at current time $t$.

1. External parameters $p$

The external parameters of MPC include the curvature of desired path $p_1=\kappa$ and desired longitudinal velocity of vehicle in the body-fixed coordinate system $p_2=v_{\mathrm{xd}}$. Since this study only focuses on the path-tracking controller design, it can be simply assumed the desired path and the curvature of desired path are already known. The desired longitudinal velocity is determined by an optimal velocity profile which can be pre-calculated offline or online as part of the motion planning algorithm. Note that the optimal velocity profile calculated offline and applied in this study can make the autonomous vehicle travel along the desired path with maximum acceleration and minimum time [32].

The MPC design has $N$ steps prediction horizon, determined as: (3) ${\kappa }_k=\kappa (X_d(t)+v_{x}k{T}_s)k=0,1,\dots N-1$(3) where $X_d(t)$ is the vehicle current time ($t$) travelling distance along the desired path. $v_x$ is the feedback vehicle velocity and $T_s$ is the sampling time of discrete-time vehicle model (2a).

1. Nonlinear prediction model

The nonlinear prediction model (2a) is determined by converting the nonlinear continuous dynamics model (1) into discrete form: (4a) $\begin{array}{rl}{v}_{x,k+1}& =v_{x,k}+T_s\left(v_{y,k}{r}_k+\frac{F_{\mathrm{xd},k}}{m}\right)\end{array}$(4a) (4b) $\begin{array}{rl}{v}_{y,k+1}& =v_{y,k}+T_s\left(-v_{x,k}{r}_k+\frac{F_{\mathrm{yd},k}}{m}\right)\end{array}$(4b) (4c) $\begin{array}{rl}{r}_{k+1}& =r_k+T_s\frac{M_{\mathrm{zd},k}}{I_z}\end{array}$(4c) (4d) $\begin{array}{rl}{F}_{\mathrm{xd},k+1}& =F_{\mathrm{xd},k}+T_s{\dot{F}}_{\mathrm{xd}}\end{array}$(4d) (4e) $\begin{array}{rl}{F}_{\mathrm{yd},k+1}& =F_{\mathrm{yd},k}+T_s{\dot{F}}_{\mathrm{yd}}\end{array}$(4e) (4f) $\begin{array}{rl}{M}_{\mathrm{zd},k+1}& =M_{\mathrm{zd},k}+T_s{\dot{M}}_{\mathrm{zd}}\end{array}$(4f) (4g) $\begin{array}{rl}{\psi }_{e,k+1}& ={\psi }_{e,k}+T_s(r_k-p_{1,k}{v}_{x,k})\end{array}$(4g) (4h) $\begin{array}{rl}{Y}_{e,k+1}& =Y_{e,k}+T_s(v_{x,k}\sin {\psi }_{e,k}+v_{y,k}\cos {\psi }_{e,k})\end{array}$(4h) where the model states $x_k=[v_{x,k};v_{y,k};r_k;F_{\mathrm{xd},k};F_{\mathrm{yd},k};M_{\mathrm{zd},k};{\psi }_{e,k};Y_{e,k}{]}^T$, and the control inputs $u_k=[{\dot{F}}_{\mathrm{xd},k};{\dot{F}}_{\mathrm{yd},k};{\dot{M}}_{\mathrm{zd},k}{]}^k$. It should be noted that in (4) the control inputs have been revised as the change rate of total desired longitudinal tire force, lateral force and yaw moment to reduce the oscillation of control inputs and improve the stability of control performance.

1. Cost function

The control targets of proposed MPC include the minimisation of yaw angle error and lateral deviation from desired path and achieving the desired target speed. The detailed mathematical presentation of the cost function can be presented as following: (5) $J=\sum _{k=0}^{N-1}{Q}_v(v_{\mathrm{xd},k}-v_{x,k}{)}^2+Q_{\psi }({\psi }_{e,k}{)}^2+Q_Y(Y_{e,k}{)}^2+{u_k}^T{Q}_u{u}_k+Q_{s1}{s}_1^2+Q_{s2}{s}_2^2$(5) where $Q_v$, $Q_{\psi }$, $Q_Y$ are scaling factors corresponding to the cost function terms of the longitudinal velocity tracking error, yaw angle tracking error and lateral deviation tracking error. $Q_u$ is a 3-by-3 diagonal matrix which include the relative penalties on the control effort. $Q_{s1}$ and $Q_{s2}$ are penalties on the soft constrains $s_1$ and $s_2$. In order to achieve good control performance, all the scaling factors in cost function are normalised and written in discrete-time form: (6a) $\begin{array}{rl}{Q}_v& =\frac{q_v{T}_s}{{(v_{\mathrm{xm}})}^2}\end{array}$(6a) (6b) $\begin{array}{rl}{Q}_{\psi }& =\frac{q_{\psi }{T}_s}{{({\psi }_{\mathrm{em}})}^2}\end{array}$(6b) (6c) $\begin{array}{rl}{Q}_Y& =\frac{q_Y{T}_s}{{(Y_{\mathrm{em}})}^2}\end{array}$(6c) (6d) $\begin{array}{rl}{Q}_u& =\left[\begin{array}{ccc}\frac{q_{u1}}{{(d{F}_{\mathrm{xm}})}^2}& 0& 0\\ 0& \frac{q_{u2}}{{(d{F}_{\mathrm{ym}})}^2}& 0\\ 0& 0& \frac{q_{u3}}{{(d{M}_{\mathrm{zm}})}^2}\end{array}\right]\end{array}$(6d) (6e) $\begin{array}{rl}{Q}_{s1}& =\frac{q_{s1}{T}_s}{{(s_{1m})}^2}\end{array}$(6e) (6f) $\begin{array}{rl}{Q}_{s2}& =\frac{q_{s2}{T}_s}{{(s_{2m})}^2}\end{array}$(6f) where $v_{\mathrm{xm}}$, ${\psi }_{\mathrm{em}}$ and $Y_{\mathrm{em}}$ are nominal values of longitudinal velocity, yaw angle error and lateral deviation, respectively. $d{F}_{\mathrm{xm}},d{F}_{\mathrm{ym}},d{M}_{\mathrm{zm}}$ are nominal values of change rate of desired longitudinal force, desired lateral force and yaw moment. $s_{1m}$ and $s_{2m}$ are nominal values of soft constraints. $q_v,q_{\psi },q_Y,q_{u1-u3},q_{s1-s2}$ are scaling values of each term, which need to be tuned as the following procedure. First the weighting of penalties terms $q_{u1-u3}$, $q_{s1-s2}$ are tuned and fixed. Then $q_v,q_{\psi },q_Y$ related to velocity error, yaw angle error and lateral deviation can be tuned individually based on the relative importance of individual term.

1. Constraints

Equations (2b) suggests that states should satisfy certain constraints. Specifically, the total desired longitudinal and lateral tire force should satisfy the friction circle limit: (7) $F_{\mathrm{xd},k}^2+F_{\mathrm{yd},k}^2-s_1^2\le (\mathrm{\mu mg}{)}^2$(7) where $\mu$ is the tire-road friction coefficient and $g$ is the constant of gravitational acceleration. The tyre-road friction coefficient can be assumed to be known through friction coefficient estimators [33].

The desired yaw moment is also constrained by its upper and lower limit: (8a) $\begin{array}{rl}& M_{\mathrm{zd},k}-M_{\mathrm{zmax}}-s_2^2\le 0\end{array}$(8a) (8b) $\begin{array}{rl}& -M_{\mathrm{zd},k}-M_{\mathrm{zmin}}-s_2^2\le 0\end{array}$(8b) The total desired yaw moment includes the yaw moment induced by the longitudinal force $M_{\mathrm{xz}}$ and lateral force $M_{\mathrm{yz}}$. Since longitudinal and lateral vehicle forces are coupled by friction limit, it is quite hard to determine the analytical solution of total yaw moment limit. Thus, in this study, a look-up table is built to explore the whole working envelope and determine the yaw moment limit, which is presented in Table 1. In Table 1, ${\delta }_{\mathrm{sf}}$ and ${\delta }_{\mathrm{sr}}$ are the searching increments of front-wheel and rear-wheel steering angle. It is also noticed the soft constraints $s_1$ and $s_2$ are applied in constraints (7, 8) to prevent the optimisation solver reaching an infeasible solution.

Table 1. Algorithm to calculate the yaw moment limit.

Assume feedback values $v_x,v_y,r,a_x,a_y,\mu ,F_{\mathrm{zi}}(\text{For i}=\text{fl},\text{fr},\text{rl},\text{rr})$ are already available

1.  For ${\delta }_f(i)=[{\delta }_{\mathrm{fmin}}:{\delta }_{\mathrm{sf}}:{\delta }_{\mathrm{fmax}}]$

2.   For ${\delta }_r(j)=[{\delta }_{\mathrm{rmin}}:{\delta }_{\mathrm{sr}}:{\delta }_{\mathrm{rmax}}]$

3.    If $v_x=0$

4.     $F_{\mathrm{yfl}}=F_{\mathrm{yfr}}=F_{\mathrm{yrl}}=F_{\mathrm{yrr}}=0$

5.    Else

6.      Calculate the side-slip angle of individual wheel $s_{\mathrm{yi}}$

7.      Calculate the lateral tire force of individual wheel $F_{\mathrm{yi}}$

8.     End

9.     $M_{\mathrm{yz}}=l_f(F_{\mathrm{yfl}}+F_{\mathrm{yfr}})-l_r(F_{\mathrm{yrl}}+F_{\mathrm{yrr}})$

10.    Calculate the maximum available individual longitudinal force

$F_{\mathrm{xi}}^{\max }=\sqrt{{(\mu F_{\mathrm{zi}})}^2-{F_{\mathrm{yi}}}^2}$

11.     Calculate $M_{\mathrm{xz}}^{\min }$ and $M_{\mathrm{xz}}^{\max }$ according to $F_{\mathrm{xi}}^{\max }$

12.     $M_z^{\min }(i)(j)=M_{\mathrm{xz}}^{\min }+M_{\mathrm{yz}},M_z^{\max }(i)(j)=M_{\mathrm{xz}}^{\max }+M_{\mathrm{yz}}$

13.    End

14.  End

15.  $M_z^{\min }=\min [M_z^{\min }(i)(j)]$, $M_z^{\max }=\max [M_z^{\max }(i)(j)]$

In order to improve the computational efficiency and satisfy the needs of real-time optimisation, the above ΝMPC problem can be linearised and solved using the fast NLP solver developed by Forces Pro [34], the detailed solving method is explained in Appendix B.

2.3 The high-level feedback path-tracking controller

In order to have the comparative study with our proposed MPC strategy, a high-level feedback path-tracking controller is presented here to calculate the desired total tire force and yaw moment [35]. The desired total longitudinal tire force can be calculated as: (9a) $F_{\mathrm{xd}}=m(-r{v}_y+{\dot{v}}_{\mathrm{xd}}-k_1{e}_1)$(9a) where $e_1=v_x-v_{\mathrm{xd}}$, $k_1>0$. ${\dot{v}}_{\mathrm{xd}}$ is the desired longitudinal acceleration and can be assumed as zero. Since the specific longitudinal acceleration target is not set in our study, we just set the desired longitudinal acceleration as zero to minimise the velocity change and improve the comfort during the desired speed tracking. The desired total lateral tire force can be calculated as: (9b) $\begin{array}{rl}{F}_{\mathrm{yd}}& =\frac{m}{\cos {\psi }_e}[-({\dot{v}}_x\sin {\psi }_e+(r-{\kappa }_d{v}_x)(v_x\cos {\psi }_e-v_y\sin {\psi }_e))\\ & +v_x\mathrm{rcos}{\psi }_e-k_2(v_x\sin {\psi }_e+v_y\cos {\psi }_e)-k_3{Y}_e]\end{array}$(9b) where $k_2>0$ and $k_3>0$. The desired yaw moment can be calculated as: (9c) $M_{\mathrm{zd}}=[-k_4(r-{\kappa }_d{v}_x)-k_5{\psi }_e+{\dot{r}}_d]I_z$(9c) where ${\dot{r}}_d$ is the desired yaw acceleration, which is assumed as 0 to improve the yaw stability of the vehicle. $k_4>0$ and $k_5>0$.

3. Low-level control allocation design

After the virtual control input $u=[F_{\mathrm{xd}};F_{\mathrm{yd}};M_{\mathrm{zd}}]$ has been determined by the upper-level MPC, the low-level CA optimisation algorithm is designed to control and allocate the individual actuator to achieve the desired virtual control force. In this study, the low-level CA algorithm is designed assuming a steer-by-wire and brake-by-wire architecture. The powertrain system includes a motor with an open differential on the front axle and two hub motors attached to the rear wheels. Two independent hydraulic brake channels are also present on the front and rear axle, respectively. Finally, the steering system can control the front- and rear-wheel steering independently.

3.1 Vehicle dynamics model for low-level CA

In this study, the design of low-level CA method are based on a two-track vehicle dynamics model. This model is more complex than the rigid-body dynamics model (1) and considers the individual wheel tire forces [36]. Note that the suspension and roll and pitch dynamics are neglected, a reasonable assumption for the on-road autonomous vehicle used in this paper.

The longitudinal motion is defined by: (10a) $m{\dot{v}}_x=m{v}_{y}r+\sum _{\begin{array}{c}i=f,r\\ j=l,r\end{array}}{F}_{\mathrm{xij}}$(10a) the lateral motion: (10b) $m{\dot{v}}_y=-m{v}_{x}r+\sum _{\begin{array}{c}i=f,r\\ j=l,r\end{array}}{F}_{\mathrm{yij}}$(10b) and the yaw motion: (10c) $I_z\dot{r}=l_f(F_{\mathrm{yfl}}+F_{\mathrm{yfr}})-l_r(F_{\mathrm{yrl}}+F_{\mathrm{yrr}})+\frac{b_f}{2}(F_{\mathrm{xfl}}-F_{\mathrm{xfr}})+\frac{b_r}{2}(F_{\mathrm{xrl}}-F_{\mathrm{xrr}})$(10c) where $F_{\mathrm{xij}}$ and $F_{\mathrm{yij}}$ are longitudinal tire force and lateral tire force of individual wheel in the body-fixed coordinate system, $\mathrm{ij}=\mathrm{fl},\mathrm{fr},\mathrm{rl},\mathrm{rr}$ represent the front left, front right, rear left and rear right wheels respectively, $l_f$ and $l_r$ are the front and rear wheel base lengths, while $b_f$ and $b_r$ are the front and rear track widths.

The longitudinal and lateral tire force of individual wheel can be determined by the following equations: (11a) $\begin{array}{rl}{F}_{\mathrm{xij}}& =F_{\mathrm{tij}}\cos {\delta }_{\mathrm{ij}}-F_{\mathrm{sij}}\sin {\delta }_{\mathrm{ij}}\end{array}$(11a) (11b) $\begin{array}{rl}{F}_{\mathrm{yij}}& =F_{\mathrm{tij}}\sin {\delta }_{\mathrm{ij}}+F_{\mathrm{sij}}\cos {\delta }_{\mathrm{ij}}\end{array}$(11b) where $F_{\mathrm{tij}}$ and $F_{\mathrm{sij}}$ are tire traction/brake force and side force of individual wheel in the wheel-centred coordinate system. For the front and rear wheel independent steering vehicle, ${\delta }_{\mathrm{fl}}={\delta }_{\mathrm{fr}}={\delta }_f$ and ${\delta }_{\mathrm{rl}}={\delta }_{\mathrm{rr}}={\delta }_r$. The individual lateral tire force can be calculated by the following equations: (12a) $\begin{array}{rl}{F}_{\mathrm{sfl}}& =C_{\mathrm{\alpha f}}{\alpha }_{\mathrm{fl}}=C_{\mathrm{\alpha f}}\left[{\delta }_f-{\tan }^{-1}\left(\frac{l_{f}r+v_y}{v_x-0.5b_{f}r}\right)\right]\end{array}$(12a) (12b) $\begin{array}{rl}{F}_{\mathrm{sfr}}& =C_{\mathrm{\alpha f}}{\alpha }_{\mathrm{fr}}=C_{\mathrm{\alpha f}}\left[{\delta }_f-{\tan }^{-1}\left(\frac{l_{f}r+v_y}{v_x+0.5b_{f}r}\right)\right]\end{array}$(12b) (12c) $\begin{array}{rl}{F}_{\mathrm{srl}}& =C_{\mathrm{\alpha r}}{\alpha }_{\mathrm{rl}}=C_{\mathrm{\alpha r}}\left[{\delta }_r+{\tan }^{-1}\left(\frac{l_{r}r-v_y}{v_x-0.5b_{r}r}\right)\right]\end{array}$(12c) (12d) $\begin{array}{rl}{F}_{\mathrm{srr}}& =C_{\mathrm{\alpha r}}{\alpha }_{\mathrm{rr}}=C_{\mathrm{\alpha r}}\left[{\delta }_r+{\tan }^{-1}\left(\frac{l_{r}r-v_y}{v_x+0.5b_{r}r}\right)\right]\end{array}$(12d) where $C_{\mathrm{\alpha f}}$ and $C_{\mathrm{\alpha r}}$ are cornering stiffness of front tire and rear tire, respectively. The individual longitudinal tire force can be determined by the wheel dynamics equation, assuming the wheel acceleration is zero: (13) $F_{\mathrm{tij}}=\frac{T_{\mathrm{ij}}}{R_{\omega }}$(13) where $R_{\omega }$ is the wheel radius and $T_{\mathrm{ij}}$ is the traction or brake torque of each wheel.

3.2 Low-level CA design

The proposed low-level CA includes two layers: in the first layer, the individual longitudinal force and steering angle are optimally distributed according to desired virtual control input; in the second layer, the distributed individual longitudinal tire forces are mapped to individual motor torque and hydraulic brake torque. When the vehicle is braking, in terms of the energy saving, we would like to try to use the negative motor torque as much as possible. After the negative motor torque reaches the maximum value, the friction brake is used to compensate the remaining required negative brake torque. This mapping strategy in the second layer can be easily achieved by a simple rule-based method. We designed this hierarchical structure to separate the first layer QCQP optimisation and the second layer rule-based method to improve the computational efficiency.

1. The first layer of CA method

According to the total virtual control force and based on the vehicle dynamics model (10)-(13), the CA equation can be presented as followings: (14a) $\begin{array}{rl}& \sum _{i=1,2,3}{u}_i=F_{\mathrm{xd}}\end{array}$(14a) (14b) $\begin{array}{rl}& C_{\mathrm{\alpha f}}(2u_4-{\alpha }_{\mathrm{fld}}-{\alpha }_{\mathrm{frd}})+C_{\mathrm{\alpha r}}(2u_5+{\alpha }_{\mathrm{rld}}+{\alpha }_{\mathrm{rrd}})=F_{\mathrm{yd}}\end{array}$(14b) (14c) $\begin{array}{rl}& -\frac{1}{2}{b}_r{u}_2+\frac{1}{2}{b}_r{u}_3+l_f{C}_{\mathrm{\alpha f}}(2u_4-{\alpha }_{\mathrm{fld}}-{\alpha }_{\mathrm{frd}})-l_r{C}_{\mathrm{\alpha r}}(2u_5+{\alpha }_{\mathrm{rld}}+{\alpha }_{\mathrm{rrd}})=M_{\mathrm{zd}}\end{array}$(14c) where the allocated control actuator input $u_1=F_{\mathrm{tf}}$, $u_2=F_{\mathrm{trl}}$, $u_3=F_{\mathrm{trr}}$, $u_4={\delta }_f$ and $u_5={\delta }_r$. It can be noticed from Equation (13(4d) $\begin{array}{rl}{F}_{\mathrm{xd},k+1}& =F_{\mathrm{xd},k}+T_s{\dot{F}}_{\mathrm{xd}}\end{array}$(4d) ) that the steering angle in tire force mapping Equation (10(4a) $\begin{array}{rl}{v}_{x,k+1}& =v_{x,k}+T_s\left(v_{y,k}{r}_k+\frac{F_{\mathrm{xd},k}}{m}\right)\end{array}$(4a) ) is considered as very small value, so $F_{\mathrm{xi}}\approx F_{\mathrm{ti}}$ and $F_{\mathrm{yi}}\approx F_{\mathrm{si}}$. ${\alpha }_{\mathrm{fld}}$,${\alpha }_{\mathrm{frd}}$,${\alpha }_{\mathrm{rld}}$,${\alpha }_{\mathrm{rrd}}$ can be determined from (11): (15) ${\alpha }_{\mathrm{id}}=f({\delta }_f,{\delta }_r,v_{\mathrm{xd}},v_{\mathrm{yd}},r_d)$(15) where $v_{\mathrm{yd}}$ is the desired lateral velocity and is assumed zero to minimise the undesired side-slip angle and $v_{\mathrm{xd}}$ and $r_d$ are the desired longitudinal and yaw rate, as determined in the high-level MPC.

In addition, the control actuation input $u_{1-5}$ should also satisfy the actuation constraints. First, the individual control actuator has its physical constraint, which is mathematically presented as follows: (16a) $\begin{array}{rl}& u_{1-3}\le \frac{T_{\mathrm{dmax}}}{R_{\omega }}\end{array}$(16a) (16b) $\begin{array}{rl}& {\theta }_{\min }\le u_4\le {\theta }_{\max }\end{array}$(16b) (16c) $\begin{array}{rl}& {\theta }_{\min }\le u_5\le {\theta }_{\max }\end{array}$(16c) where $T_{\mathrm{dmax}}$ represents the maximum driving toque of the individual motor, and ${\theta }_{\min }$, ${\theta }_{\max }$ are the minimum and maximum steering angle respectively. Then the friction circle constraint should also be considered to show the maximum available friction force of each tire: (17a) $\begin{array}{rl}& 0.25u_1^2+C_{\mathrm{\alpha f}}^2(u_4-{\alpha }_{\mathrm{fld}}{)}^2\le (\mu F_{\mathrm{zfl}}{)}^2\end{array}$(17a) (17b) $\begin{array}{rl}& 0.25u_1^2+C_{\mathrm{\alpha f}}^2(u_4-{\alpha }_{\mathrm{frd}}{)}^2\le (\mu F_{\mathrm{zfr}}{)}^2\end{array}$(17b) (17c) $\begin{array}{rl}& u_2^2+C_{\mathrm{\alpha r}}^2(u_5+{\alpha }_{\mathrm{rld}}{)}^2\le (\mu F_{\mathrm{zrl}}{)}^2\end{array}$(17c) (17d) $\begin{array}{rl}& u_3^2+C_{\mathrm{\alpha r}}^2(u_5+{\alpha }_{\mathrm{rrd}}{)}^2\le (\mu F_{\mathrm{zrr}}{)}^2\end{array}$(17d) where $F_{\mathrm{zfl}}$, $F_{\mathrm{zfr}}$, $F_{\mathrm{zrl}}$, $F_{\mathrm{zrr}}$ are the vertical loads on the wheels. Considering load transfer effects, the vertical load of each wheel can be estimated by the following equations [37]: (18a) $\begin{array}{rl}{F}_{\mathrm{zfl}}& =\frac{m}{l_f+l_r}\left(\frac{1}{2}g{l}_r-\frac{1}{2}({\dot{v}}_x-v_{y}r)h-\frac{l_r}{b_f}({\dot{v}}_y+v_{x}r)h\right)\end{array}$(18a) (18b) $\begin{array}{rl}{F}_{\mathrm{zfr}}& =\frac{m}{l_f+l_r}\left(\frac{1}{2}g{l}_r-\frac{1}{2}({\dot{v}}_x-v_{y}r)h+\frac{l_r}{b_f}({\dot{v}}_y+v_{x}r)h\right)\end{array}$(18b) (18c) $\begin{array}{rl}{F}_{\mathrm{zrl}}& =\frac{m}{l_f+l_r}\left(\frac{1}{2}g{l}_f+\frac{1}{2}({\dot{v}}_x-v_{y}r)h-\frac{l_f}{b_r}({\dot{v}}_y+v_{x}r)h\right)\end{array}$(18c) (18d) $\begin{array}{rl}{F}_{\mathrm{zrr}}& =\frac{m}{l_f+l_r}\left(\frac{1}{2}g{l}_f+\frac{1}{2}({\dot{v}}_x-v_{y}r)h+\frac{l_f}{b_r}({\dot{v}}_y+v_{x}r)h\right)\end{array}$(18d) where $h$ is the distance of the vehicle Centre of Gravity (CG) from the ground. It is noted that the accelerations ${\dot{v}}_x$ and ${\dot{v}}_y$ and yaw rate $r$ can be obtained from an Inertial Measurement Unit (IMU), and the measurement signals of $v_x$ and $v_y$ are available through the Global Positioning System (GPS) unit.

In order to solve the CA problem (14-18), the optimisation problem can be designed to minimise the virtual input error and the cost of using the control actuators: (19a) $\underset{u_i(\text{for }i=1,2,3,4,5)}{min}J=({\tau }_{\mathrm{RHS}}-\mathrm{Bu}(k+j){)}^T{Q}_1({\tau }_{\mathrm{RHS}}-\mathrm{Bu}(k+j))+u(k+j{)}^T{Q}_{2}u(k+j)$(19a) The cost function (19a) can be re-written in the following QCQP form: (19b) $\underset{u_i(\text{for }i=1,2,3,4,5)}{min}J\approx u(k+j{)}^T(B^T{Q}_{1}B+Q_2)u(k+j)+(-2B{Q}_1{\tau }_{\mathrm{RHS}})u(k+j)$(19b) s.t. (19c) $\begin{array}{rl}& {\left[\begin{array}{c}{u}_1\\ u_4\end{array}\right]}^T\left[\begin{array}{cc}0.25& 0\\ 0& C_{\mathrm{\alpha f}}^2\end{array}\right]\left[\begin{array}{c}{u}_1\\ u_4\end{array}\right]+\left[\begin{array}{cc}0& -2{\alpha }_{\mathrm{fld}}{C}_{\mathrm{\alpha f}}^2\end{array}\right]\left[\begin{array}{c}{u}_1\\ u_4\end{array}\right]\le (\mu F_{\mathrm{zfl}}{)}^2-(C_{\mathrm{\alpha f}}{\alpha }_{\mathrm{fld}}{)}^2\end{array}$(19c) (19d) $\begin{array}{rl}& {\left[\begin{array}{c}{u}_1\\ u_4\end{array}\right]}^T\left[\begin{array}{cc}0.25& 0\\ 0& C_{\mathrm{\alpha f}}^2\end{array}\right]\left[\begin{array}{c}{u}_1\\ u_4\end{array}\right]+\left[\begin{array}{cc}0& -2{\alpha }_{\mathrm{frd}}{C}_{\mathrm{\alpha f}}^2\end{array}\right]\left[\begin{array}{c}{u}_1\\ u_4\end{array}\right]\le (\mu F_{\mathrm{zfr}}{)}^2-(C_{\mathrm{\alpha f}}{\alpha }_{\mathrm{frd}}{)}^2\end{array}$(19d) (19e) $\begin{array}{rl}& {\left[\begin{array}{c}{u}_2\\ u_5\end{array}\right]}^T\left[\begin{array}{cc}1& 0\\ 0& C_{\mathrm{\alpha r}}^2\end{array}\right]\left[\begin{array}{c}{u}_2\\ u_5\end{array}\right]+\left[\begin{array}{cc}0& 2{\alpha }_{\mathrm{rld}}{C}_{\mathrm{\alpha r}}^2\end{array}\right]\left[\begin{array}{c}{u}_2\\ u_5\end{array}\right]\le (\mu F_{\mathrm{zrl}}{)}^2-(C_{\mathrm{\alpha r}}{\alpha }_{\mathrm{rld}}{)}^2\end{array}$(19e) (19f) $\begin{array}{rl}& {\left[\begin{array}{c}{u}_3\\ u_5\end{array}\right]}^T\left[\begin{array}{cc}1& 0\\ 0& C_{\mathrm{\alpha r}}^2\end{array}\right]\left[\begin{array}{c}{u}_3\\ u_5\end{array}\right]+\left[\begin{array}{cc}0& 2{\alpha }_{\mathrm{rrd}}{C}_{\mathrm{\alpha r}}^2\end{array}\right]\left[\begin{array}{c}{u}_3\\ u_5\end{array}\right]\le (\mu F_{\mathrm{zrr}}{)}^2-(C_{\mathrm{\alpha r}}{\alpha }_{\mathrm{rrd}}{)}^2\end{array}$(19f) where $u^T=\left[\begin{array}{ccc}\begin{array}{cc}{u}_1& u_2\end{array}& \begin{array}{cc}{u}_3& u_4\end{array}& u_5\end{array}\right],B=\left[\begin{array}{cc}\begin{array}{cc}1& 1\end{array}& \begin{array}{ccc}1& 0& 0\end{array}\\ \begin{array}{cc}0& 0\end{array}& \begin{array}{ccc}0& 2C_{\mathrm{\alpha f}}& 2C_{\mathrm{\alpha r}}\end{array}\\ \begin{array}{cc}0& -0.5b_r\end{array}& \begin{array}{ccc}0.5b_r& 2l_f{C}_{\mathrm{\alpha f}}& 2l_r{C}_{\mathrm{\alpha r}}\end{array}\end{array}\right]$ ${\tau }_{\mathrm{RHS}}=\left[\begin{array}{c}{F}_x\\ F_y+C_{\mathrm{\alpha f}}({\alpha }_{\mathrm{fld}}+{\alpha }_{\mathrm{frd}})-C_{\mathrm{\alpha r}}({\alpha }_{\mathrm{rld}}+{\alpha }_{\mathrm{rrd}})\\ M_z+l_f{C}_{\mathrm{\alpha f}}({\alpha }_{\mathrm{fld}}+{\alpha }_{\mathrm{frd}})+l_r{C}_{\mathrm{\alpha r}}({\alpha }_{\mathrm{rld}}+{\alpha }_{\mathrm{rrd}})\end{array}\right],$

Note also that the boundary constraint (16) should be satisfied. $Q_1$ and $Q_2$ are the diagonal weighting matrices on the virtual control error and control effort respectively. $Q_1=\text{diag}[q_{11};q_{12};q_{13}]$, $q_{11}>0,q_{12}>0,q_{13}>0$, where $q_{11}$ is the weighting factor for achieving the desired total longitudinal tire force, $q_{12}$ is the weighting factor for achieving the desired total lateral tire force, $q_{13}$ is the weighting factor for achieving the desired yaw moment; $Q_2=\text{diag}[q_{21};q_{22};q_{23};q_{24};q_{25}]$, $q_{21}>0$, $q_{22}>0$, $q_{23}>0$, $q_{24}>0$, $q_{25}>0$, where $q_{21}$ is the weighting factor for penalising the allocated total front tire longitudinal force, $q_{22}$ is the weighting factor for penalising the allocated rear left tire longitudinal force, $q_{23}$ is the weighting factor for penalising the allocated rear right tire longitudinal force, $q_{24}$ is the weighting factor for penalising the allocate front wheel steering angle, $q_{25}$ is the weighting factor for penalising the allocate rear wheel steering angle. In order to achieve the good control performance, the weightings of $q_{11},q_{12},q_{13}$ are set as much bigger than the weightings of $q_{21},q_{22},q_{23},q_{24},q_{25}$.

In order to improve the computational efficiency and satisfy the needs for real-time CA optimisation, the above QCQP problem can be solved by the time-efficient solver based on Primal Dual Interior Point (PDIP) method in Forces Pro [38] in appendix A.

1. The second layer of CA method

In the second layer of CA, the allocated individual longitudinal force from first layer is mapped into the individual motor torque or hydraulic brake torque. The mapping strategies of front axle and rear axle are presented in Figure 3(a,b), respectively. The main idea of the mapping strategy is trying to use the electric motor as much as possible (as electric motor can generate both traction torque and regenerated brake torque). If the electric motor cannot generate enough brake torque, the hydraulic brake will compensate the brake torque. In Figure 3(a), $T_{\mathrm{mf}}$ and $T_{\mathrm{bf}}$ represent the electric motor torque and hydraulic brake torque on the front axle respectively. $T_{\mathrm{mf}}^{\min }$ and $T_{\mathrm{mf}}^{\max }$ are the motor torque limit. In Figure 3(b), $T_{\mathrm{mrl}}$ and $T_{\mathrm{mrr}}$ are electric motors of rear left and rear right wheel and $T_{\mathrm{br}}$ is the hydraulic brake torque on the rear axle. $T_{\mathrm{mrl}}^{\min }$,$T_{\mathrm{mrl}}^{\max }$, $T_{\mathrm{mrr}}^{\min }$ and $T_{\mathrm{mrr}}^{\max }$ are the motor torque limits. It is noted that constraints (16a) in first layer of CA has guaranteed $u_1\le \frac{T_{\mathrm{mf}}^{\max }}{R_{\omega }}$, $u_2\le \frac{T_{\mathrm{mrl}}^{\max }}{R_{\omega }}$ and $u_3\le \frac{T_{\mathrm{mrr}}^{\max }}{R_{\omega }}$, which is the precondition in Figure 3.

Figure 3. The allocation of motor torque and hydraulic brake toque (a) front axle (b) rear axle.

4. Simulation results

In this section, the proposed hierarchical path-tracking controller is validated using a high-fidelity vehicle dynamics model in CarMaker. The vehicle model is built based on an autonomous test platform developed by Delta Motorsport Limited. Specifically, the software simulation environment to test the optimisation control strategy is the jointed simulation of Matlab Simulink, Forces Pro optimisation solver and IPG CarMaker. The vehicle plant model and simulation environment (road and surrounding environment) are implemented in IPG CarMaker, and the high-level MPC controller and low-level control allocation method are implemented based on Matlab Simulink and Forces Pro optimisation solver (a Simulink optimisation code package). The IPG CarMaker based vehicle plant model and Simulink based optimisation strategy are integrated together through the CarMaker-Simulink co-simulation. The computer hardware to run the simulation is a IBM Thinkpad laptop (Intel(R) Core(TM) i7-5500U CPU @ 2.40 GHz and 8 GB RAM).

The simulation scenario chosen on the Silverstone racetrack, and the autonomous vehicle is trying to track the road center-line. The velocity profile with minimum-time solution along the whole racetrack is pre-calculated offline. Both the desired path and velocity profile are inputted into the hierarchical controller, and the controller distributes the steering and driving control inputs into the CarMaker vehicle model. The CarMaker model also provides the necessary information on the state values which can be then fed back to the hierarchical controller. The main vehicle parameters are presented in Table 2.

Table 2. Main parameters of Carmaker model.

Symbol	Meaning	Values
$m$	Vehicle total mass	700.28 kg
$I_z$	Moment inertial around yaw axle	1597.717 kg.m^2
$l_f$	Front wheel base	0.999 m
$l_r$	Rear wheel base	0.996 m
$b_f$	Front track width	1.52 m
$b_r$	Rear track width	1.52 m
$R_{\omega }$	Wheel radius	0.32 m
$C_{\alpha }$	Cornering stiffness for front and rear tyres	29220 N/rad
$\mu$	Tire-road friction coefficient	1
$k_{1-5}$	Control gains of high-level feedback controller	[3 20 5 20 200]

A particular s-curve section of Silverstone track is selected to present the control performance of proposed controller. In order to find out the best sampling time and prediction horizon of the proposed high-level path-tracking MPC, the sensitivity analysis on different sampling time and prediction horizons (sampling time = 0.05, 0.075 s; prediction time = 2s, 4s, 6s) is carried out by comparing the path-tracking performance and velocity-tracking performance. When sampling time and prediction time are 0.05 and 2 s. respective, the best control performance can be achieved. Thus, the simulation results below are presented by using sampling and prediction time of 0.05 and 2 s. The control performance of the proposed CA (‘MPC + CA’) in the low-level is compared with two other, simpler vehicle configurations: the CA method without rear-wheel steering is simply presented as ‘MPC + CA without RS’ and the CA method without torque vectoring is presented as ‘MPC + CA without TV’. It is noted that for the ‘MPC + CA without RS’, the weighting value related to rear-wheel steering angle in $Q_2$ is adjusted to disable the rear-wheel steering in (19). For the ‘MPC + CA without TV’, the B matrix is adjusted to assume the motor torque is equally distributed among two rear wheels in (19). In order to show the advantage of proposed high-level MPC, the simulation results of a high-level feedback path-tracking controller in (8) and proposed low-level CA method are also presented as ‘FB + CA’ in following figures. Figure 4(a) shows the desired path. Figure 4(b) presents the path-tracking performance and Figure 4(c) shows velocity tracking performance of different methods. According to Figure 4(b), ‘MPC + CA’ has the smallest lateral tracking error compared with ‘MPC + CA without TV’, ‘MPC + CA without RS’ and ‘FB + CA’. Figure 4(c) demonstrates that ‘MPC + CA’ and ‘FB + CA’ can fast track the desired longitudinal velocity profile though the online optimisation and achieve the maximum cornering speed, while other two methods both show a certain level of delay. Figure 5 shows the total normalised acceleration, which is calculated by the following equation: (20) $\overline{a}=\frac{\sqrt{a_x^2+a_y^2}}{\mathrm{\mu g}}$(20) where $a_x$ and $a_y$ are actual longitudinal and lateral acceleration, respectively. It can be observed from Figure 5 that the ‘MPC + CA’ and the ‘MPC + CA without TV’ have larger normalised acceleration than ‘MPC + CA without RS’, which means that CA and CA without TV can take greater usage of tire force to achieve the limit-handling performance. The normalised acceleration of ‘FB + CA’ shows big spike, which causes the sudden increase of tracking error in Figure 4(b) and suggests the less stability of high-level feedback controller compared with MPC. It should be noted that the normalised acceleration limit is larger than 1 in this study because a more ‘potent’ tire model is assumed.

Figure 4. The tracking performance in the simulation (left) desired path (middle) lateral path-tracking error (right) velocity-tracking performance.

Figure 5. The normalised acceleration in the simulation.

In Figure 6, the simulation results of robustness analysis on proposed ‘MPC + CA’ controller are presented. First, it is assumed that the friction coefficient can not be estimated accurately. Figure 6 compares the path-tracking error and longitudinal velocity tracking error when the estimated friction coefficient is 1 and 0.7, respectively (the actual value is 1). When the friction coefficient is wrongly estimated as 0.7, the path-tracking error is relatively larger. This is because the actual total acceleration is larger than the wrongly estimated friction limit in low-level CA and the control allocation performance is compromised. Nevertheless, the path-tracking error is still acceptable as it is smaller than 0.5 m. Secondly, the measurement white noise of the longitudinal acceleration (variance is 0.03 m/s²), lateral acceleration (variance is 0.06 m/s²) and yaw rate (variance is 0.0003 rad/s) are considered in the simulation and the control performance is similar as the results when the measurement noise is not considered. In general, Figure 6 suggests that the proposed ‘MPC + CA’ controller has good robustness performance.

Figure 6. The robustness analysis of the proposed control strategy

Figure 7 shows the tire forces and yaw moment performance for different controllers. It can be seen in Figure 7 that the desired longitudinal force and yaw moment are tracked quite well. The desired lateral force is much smaller than actual values. This is because the desired lateral tire force is calculated based on the lateral position tracking error, the small desired lateral force suggests that the vehicle has good lateral tracking performance. Figure 7(a) suggests that the total longitudinal tire forces of ‘MPC + CA without RS’ are much smaller than ‘MPC + CA’ and ‘MPC + CA without TV’, which can be the reason why ‘MPC + CA without RS’ have worse longitudinal velocity tracking performance compared to the other two strategies. According to Figure 7(b), the ‘MPC + CA without RS’ shows smaller total lateral tire force since the vehicle can only achieve front-wheel steering and the total available lateral tire force that can be generated is reduced. Figure 7 (d) shows that ‘MPC + CA without TV’ cannot generate additional torque-vectoring yaw moment. Figures 7(a,b,d) demonstrate that the total forces and the additional torque-vectoring yaw moment in the case of the ‘MPC + CA’ are larger than the other methods and consequently the overall path-tracking and velocity-tracking performance is the best amongst the controllers. The CA without TV disables the torque vectoring function and cannot generate additional controlled yaw moment when the vehicle is making the sharp turning as shown in Figure 7(d), so the tracking performance is compromised as shown in Figure 4. The CA without RS disables the rear-wheel steering function and limits the available total lateral tire force so the tracking performance is also compromised (as shown also in Figure 4). Finally, Figure 7(a-d) suggest that tire forces and yaw moment of ‘FB + CA’ method shows big spike during the sharp turning, which cause the compromised path-tracking performance.

Figure 7. The actual tire force and yaw moment in the simulation (a) total longitudinal force (b) total lateral force (c) total yaw moment (d) torque vectoring yaw moment.

Figure 8 compares the torque outputs of the three motors. The torques from two rear hub motors of ‘MPC + CA’ and ‘FB + CA’ in Figures 8(a,d) and ‘MPC + CA without RS’ in Figure 8(c) show obvious torque vectoring performance. Figure 9 suggests that the steering angles from the ‘MPC + CA’ and ‘MPC + CA without TV’ are larger than ‘MPC + CA without RS’, which is consistent with the results for the total lateral force in Figure 7(b). In Figure 9(c), the rear wheel steering angle of ‘MPC + CA without RS’ is zero. Figure 9(d) shows the steering angles of ‘FB + CA’ is much bigger than other methods, which impairs the stability of the vehicle.

Figure 8. The motor torque in the simulation (a) MPC + CA (b) MPC + CA without TV (c) MPC + CA without RS (d) FB + CA.

Figure 9. The wheel steering angle in the simulation (a) MPC + CA (b) MPC + CA without TV (c) MPC + CA without RS (d) FB + MPC.

In order to validate the real-time performance of the proposed hierarchical controller, the solving time of high-level path-tracking MPC and low-level CA method are presented in Figure 10. Since the sampling time of the proposed controller is 0.05s, the solving time in Figure 10 is smaller than 0.05 s and demonstrates that the real-time computational requirement can be satisfied. It is noted that the high-level MPC solver takes longer time to find a solution around 2, 5.5, 17, 19 and 21.5 s. The reason behind this is more iterations are required for the solver to find the optimal solution, or even the optimal solution cannot be found within the max iteration and only the sub-optimal solution is found. It is also noted that unlike the high-level MPC controller which required to solve the MPC optimisation problem in prediction steps, the low-level CA controller only need to solve the standard QCQP optimisation problem without the prediction model. So the optimisation solving time of low-level controller is much smaller than the high-level controller.

Figure 10. The solving time of CA in the simulation.

In the literature, the execution time of MPC path-tracking controller [39] based on vehicle lateral dynamics is 100 ms with the prediction steps of 20. Wang et al. suggests a classic Linear-Time-Varying (LPV) MPC based on a 3 DOF nonlinear vehicle dynamics model including the Brush tyre model and the analysis on the computational efficiency is carried out [40]. They set the sampling time as 50 ms (The execution time is close to 50 ms) with the time steps of 14. For our designed two-layer MPC, the execution time is around 20 ms with time steps of 40, which shows better computational efficiency compared with above two holistic MPC approaches.

5. Conclusion

In this study, a hierarchical two-level integrated path-tracking and CA controller is designed. The high-level path-tracking controller is composed by on a nonlinear rigid body dynamics model and low-level CA is designed using a QCQP formulation. The major findings can be summarised as follows:

1. By selecting the best sampling time of 0.05 s and prediction time of 2 s, and reduce the lateral tracking error within 0.5 m with maximum cornering speed during limit handling scenario in the suggested simulation platform, and shows better performance than a benchmark feedback path-tracking controller.

2. The proposed low-level CA method fully utilizes torque vectoring of two rear motors and front and rear-wheel independent steering, and the path-tracking and velocity-tracking performance are improved when compared with a controller without torque vectoring or a controller without rear-wheel steering.

3. The total solving time of high-level MPC and low-level CA is smaller than the sampling time for the whole strategy, which shows that the proposed method can be implemented in real-time.

In future work, the real-time performance will be tested on a real-time target machine (such as dSpace or Speedgoat) before deploying the final solution on the autonomous test platform.

Disclosure statement

No potential conflict of interest was reported by the author(s).

Correction Statement

This article has been republished with minor changes. These changes do not impact the academic content of the article.

Additional information

Funding

This work was supported by Innovate UK through AID-CAV project.

Appendices

Appendix A

The QCQP optimisation problem (19b) can be rewritten as: (A1a) $\underset{u_i(\text{for }i=1,2,3,4,5)}{min}J=l(u_i(k+j))$(A1a) The inequality constraints (19c–19f) can be rewritten as: (A1b) $g_j(u_i(k+j))\le 0$(A1b) In order to solve the MPC optimisation problem (A1) in a time-efficient manner, a primal-dual interior point method is applied. Particularly, the computation of search direction is achieved by solving the Karush-Kuhn-Tucker (KKT) system and the details are given as followings:

1. The formulation of KKT system

The KKT optimal condition for problem (A1) is given as: (A2a) $\begin{array}{rl}& h(u_i)+J(u_i{)}^T\lambda =0\end{array}$(A2a) (A2b) $\begin{array}{rl}& g_{j=1,2,3,4}(u_i)+s=0\end{array}$(A2b) (A2c) $\begin{array}{rl}& \wedge S=0\end{array}$(A2c) where $h(u_i)=\mathrm{\nabla l}(u_i(k+1))$, $g_j(u_i)=g_j(u_i(k+1))$,$J(u_i)=\text{blkdiag}[\mathrm{\nabla }{g}_j(u_i(k+1))]$. $\lambda$ is the Lagrange multiplier. $s$ is the slack variable and $S=\text{diag}(s)$, $\wedge =\text{diag}(\lambda )$.

1. Search direction

The Newton system is determined by linearising KKT system: (A3) $\left[\begin{array}{ccc}H(u_i,\lambda )& J{(u_i)}^T& 0\\ J{(u_i)}^T& 0& I\\ 0& S& \wedge \end{array}\right]\left[\begin{array}{c}\mathrm{\Delta }{u}_i\\ \mathrm{\Delta \lambda }\\ \mathrm{\Delta s}\end{array}\right]=-\left[\begin{array}{c}{r}_C\\ r_I\\ r_s\end{array}\right]$(A3) Where $H(u_i,\lambda )=\mathrm{\nabla h}(u_i)+\text{blkdiag}[{\mathrm{\nabla }}^2{g}_j(u_i){\lambda }_k]$.

1. Block elimination

Block elimination can be applied to simplify the search direction computation. It can be assumed that the elements of $\wedge$ are strictly positive, $\mathrm{\Delta }s$ can be eliminated from (A3) by using $\mathrm{\Delta }s={\wedge }^{-1}(r_s-S\mathrm{\Delta \lambda })$. The Lagrange multiplier can be rearranged as following: (A4) $\mathrm{\Delta }\lambda =S^{-1}\wedge (r_I+J(u_i)\mathrm{\Delta }{u}_i)+S^{-1}{r}_s$(A4) Thus, the search direction equations (B3) can be simplified as: (A5) $\mathrm{\Phi \Delta }{u}_i=-r_d$(A5) where $\mathrm{\Phi }=H(u_i,\lambda )+J(u_i{)}^T{S}^{-1}\wedge J(u_i)$, $r_d=r_C+J(u_i{)}^T{S}^{-1}\wedge r_I+J(u_i{)}^T{S}^{-1}{r}_s$. Thus, the search direction of (A5) can be solved in a time-efficient manner.

Appendix B

The nonlinear optimisation problem (2-8) can be rewritten as the followings: (B1a) $\underset{u_i,x_j}{min}J=\sum _{k=1}^N{l}_k(z_{i+j}(k))$(B1a) subject to: (B1b) $\begin{array}{rl}& c_k(z_{i+j}(k+1),z_{i+j}(k))=0\end{array}$(B1b) (B1c) $\begin{array}{rl}& h_k(z_{i+j}(k))\le 0\end{array}$(B1c) The nonlinear MPC problem (B1) can be rewritten as the following linearised KKT version: (B2) $\left[\begin{array}{cc}\begin{array}{cc}H(z_{i+j},y_c)& J_{\mathrm{eq}}{(z_{i+j})}^T\\ J_{\mathrm{eq}}(z_{i+j})& 0\end{array}& \begin{array}{cc}{J}_{\mathrm{ineq}}{(z_{i+j})}^T& 0\\ 0& 0\end{array}\\ \begin{array}{cc}{J}_{\mathrm{ineq}}(z_{i+j})& 0\\ 0& 0\end{array}& \begin{array}{cc}0& I\\ S& Y_D\end{array}\end{array}\right]\left[\begin{array}{c}\begin{array}{c}\mathrm{\Delta z}\\ \mathrm{\Delta }{y}_c\end{array}\\ \begin{array}{c}\mathrm{\Delta }{y}_d\\ \mathrm{\Delta s}\end{array}\end{array}\right]=-\left[\begin{array}{c}\begin{array}{c}{r}_S\\ r_C\end{array}\\ \begin{array}{c}{r}_D\\ r_N\end{array}\end{array}\right]$(B2) where $H(z_{i+j},y_c)={\mathrm{\nabla }}_{z_{i+j}}^2{l}_k(z_{i+j}(k))+\sum _{m=1}^{r-1}{y}_{c,k-1[m]}{\mathrm{\nabla }}_{z_{i+j}}^2{c}_{k-1[m]}(z_{i+j}(k),z_{i+j}(k-1))+\sum _{m=1}^r{y}_{c,k[m]}{\mathrm{\nabla }}_{z_i}^2{c}_{k[m]}(z_{i+j}(k+1),z_{i+j}(k))$. $y_c$ is the Lagrange multiplier of the equality constraints, and $y_d$ is the Lagrange multiplier of the inequality constraints. $s$ is the slack variable and $S=\text{diag}(s)$. $J_{\mathrm{ineq}}$ is the block diagonal consisting of $\mathrm{\nabla }{h}_k(z_{i+j})$, while $J_{\mathrm{eq}}$ is the block diagonal consisting of $\mathrm{\nabla }{c}_k(z_{i+j})$. $Y_D$ is the diagonal matrix with the elements of $y_d$ on its diagonal.

When dealing with the nonlinear MPC problem instead of convex QP problem, $H$ is a nonlinear function of $z_{i+j}$ and $y_c$, while the Jacobians of equality constraints $J_{\mathrm{eq}}$ and inequality constraints $J_{\mathrm{ineq}}$ are a function of $z_{i+j}$. It is important to mention that it is general difficult to ensure the KKT formulation (B2) has the desired inertial as QP problem, so additional care must be taken to solve (B2).

First, the linear system is reduced by eliminating $\mathrm{\Delta }s$ and $\mathrm{\Delta }{y}_d:$ (B3) $\left[\begin{array}{cc}\mathrm{\Phi }& J_{\mathrm{eq}}{(z_{i+j})}^T\\ J_{\mathrm{eq}}(z_{i+j})& 0\end{array}\right]\left[\begin{array}{c}\mathrm{\Delta }z\\ \mathrm{\Delta }{y}_c\end{array}\right]=-\left[\begin{array}{c}{r}_d\\ r_C\end{array}\right]$(B3) where $\mathrm{\Phi }=H(z_{i+j},y_c)+J_{\mathrm{ineq}}(z_{i+j}{)}^T{S}^{-1}{Y}_D{J}_{\mathrm{ineq}}(z_{i+j})$, $r_d=r_S+J_{\mathrm{ineq}}(z_{i+j}{)}^T{S}^{-1}{Y}_D{r}_D-J_{\mathrm{ineq}}(z_{i+j}{)}^T{S}^{-1}{r}_N$.

Secondly, a more compact and symmetric system can be obtained as following equation: (B4) $Y\mathrm{\Delta }{y}_c=\beta$(B4) where $Y=J_{\mathrm{eq}}(z_{i+j}){\mathrm{\Phi }}^{-1}{J}_{\mathrm{eq}}(z_{i+j}{)}^T$, $\beta =r_C-J_{\mathrm{eq}}(z_{i+j}){\mathrm{\Phi }}^{-1}{r}_d$.

Forces NLP provides linear solver to solve the symmetric form (B3) and normal form (B4). For the QP formulation, $\mathrm{\Phi }$ and $Y$ can be assumed as positive definite and the problem can be easily solved by the Cholesky decomposition. However, for the nonlinear MPC formulation, the Hessian of Lagrangian is potentially indefinite and $Y$ cannot be factorised with the Cholesky decomposition. A number of methods are applied to approximate the Hessian matrix and the positive definite can be guaranteed [38,41].