Path following of underactuated surface ships based on model predictive control with neural network

A model predictive control approach is proposed for path following of underactuated surface ships with input saturation, parameters uncertainties, and environmental disturbances. An Euler iterative algorithm is used to reduce the calculation amount of model predictive control. The matter of input saturation is addressed naturally and flexibly by taking advantage of model predictive control. The mathematical model group (MMG) model as the internal model improves the control accuracy. A radial basis function neural network is also applied to compensate the total unknowns including parameters uncertainties and environmental disturbances. The numerical simulation results show that the designed controller can force an underactuated ship to follow the desired path accurately in the case of input saturation and time-varying environmental disturbances including wind, current, and wave.


Introduction
Path following of surface ships has been a long-standing control problem that has attracted attention from the control community for many years. 1 As a matter of fact, the path following control possesses the significant position in some marine applications, such as seawater monitoring for the fixed points, laying and maintaining the underwater pipelines, even tracking control of the ordinary merchant ships. The underactuated nature of these problems, namely with more variables to be controlled than the number of control actuators, coupled with the matter of uncertain parameters, external disturbances as well as input constraint, renders the control problem both challenging and interesting. 2,3 The path following control problem of underactuated ships has been tackled in a number of research studies in the last few years. For example, a position prediction model was presented by Nagai and Watanabe, 4 which estimated positions by the integral in continuous time and predicted future positions by considering only constant velocities. In the study of Shen and Dai, 5 the iterative sliding mode was designed by a hyperbolic tangent function, and both neural network (NN) and reinforced learning were used to inhibit the chattering of the control input.
Methods to deal with the uncertain parameters or external disturbances could also be found in previous studies. In the study of Wang et al., 6 a robust and adaptive control algorithm was proposed to cope with the random disturbances and uncertain parameters, in which the Serret-Frenet frame and output-redefinition were introduced to simplify the controller design. In the study of Herman and Adamski, 7 the robustness was promoted by designing a non-adaptive velocity tracking algorithm for fully actuated vehicles. Li et al. 8 introduced an active disturbance rejection control (ADRC) controller with sliding mode which can solve the path following problem of uncertainties including internal dynamic and external disturbances. In the study of Lekkas and Fossen, 9 a sliding mode controller was employed to stabilize the heading angle, and the constant external disturbances could be handled. It is noted that the sliding mode and incremental feedback methods were used by Liu et al. 10 to avoid effectively the trajectory tracking with unknown disturbances such as wind or current. In the studies of Mu et al. 11 and Liang et al., 12 both surge control and yaw control were considered by designing the sliding mode trajectory tracking controller, and the unknown dynamics or disturbances were packaged for compensation via NN together with minimum learning parameter method (MLP). The ADRC was proposed by Huang and Fan 13 to deal with the path following problem, and the extended state observer (ESO) and adaptive laws were used to estimate the external disturbances and sideslip angle, respectively. In the study of Vo et al., 14 an adaptive autopilot based on backstepping methodology was designed to control ship heading with external disturbances. In the study of Zhang et al., 15 the nonlinear disturbance observer was designed to estimate disturbances, and a predictive controller was offered based on Serret-Frenet frame, but the advantage of predictive control in saturated suboptimal solutions was not utilized.
To solve the problems of both the uncertain parameters and environmental disturbances, the robust and adaptive control algorithm was also used by Tong and Li. 16,17 An adaptive takagi and sugeno (T-S) fuzzy NN control method was proposed by Dong et al. 18 to deal with the external disturbance, but the simulation model was very different from a real ship. Model uncertainties were approximated and compensated by virtue of the concise radial basis function (RBF) adaptive NN in the study of Zhang et al. 19 Some investigations to address the problem of input saturation could be found in previous research studies. Zheng et al. 20 proposed a linear model predictive control (MPC) and a nonlinear MPC to demonstrate the advantage of MPC on dealing with the matter of saturation. A novel robust MPC was presented by Li et al. 21 to solve the vessel's roll constraints. In the study of Shen and Jing, 22 both neuron adaptive and iterative sliding mode strategy were developed for the path tracking of underactuated ships. In the study of Yu et al., 23 the input saturation was considered via an auxiliary system and the MLP was provided to approximate the unknown dynamic without external disturbances.
A backstepping method was designed by Zheng and Feroskhan, 24 and the disturbance observer and auxiliary system were used to estimate external disturbances and cope with the input saturation, respectively. In the study of Qiu et al., 25 the external disturbances, dynamical uncertainties, and input saturation were settled effectively by the adaptive technology, MLP, and auxiliary dynamic system, respectively.
In the aforementioned studies, the problem of unknown dynamical or external disturbances could be addressed effectively. However, few have taken into account the input saturation, whereas the poor control performances or even actuator damages would happen in real implementation if these constraints were neglected in the controller design. 20 Some research studies have pointed out the solutions of the saturation such as introducing an iterative sliding mode theory 5,10,22 or an auxiliary system. [23][24][25] Comparing to the schemes above, the MPC methodology proposed in this article is more natural and flexible 20,26 as just the value of the constrained rudder angle is considered directly in case of solving a quadratic programming (QP) 27 problem. The designed MPC is simpler than those proposed in the literature, 6,15,20,21 because the Euler iterative algorithm is introduced to predict the future states. The MMG model is used to improve the prediction accuracy as the internal model. Moreover, with reference to the literature, 19,23,25,28 a RBF NN is used to compensate the uncertainties for the improvement of the internal model accuracy and the controller robustness.
The remainder of this article is organized as follows. The kinematic and dynamic model of the surface ship is presented in the next section. The MPC with RBF NN method is developed in the third section. The fourth section provides the stability analysis. Simulation results are presented and analyzed in the fifth section, followed by conclusions.

An underactuated surface ship model
The ship position in the horizontal plane and the motion parameters with ocean current disturbance are shown in Figure 1. The x-axis and y-axis point to the true North and East, respectively. ' is the heading angle from the x-axis. u and v are the surge and sway velocity over the ground, and V ¼ (u 2 þ v 2 ) 1/2 is the ship speed over ground. r is the yaw rate. V c and ' c are the ocean current speed and direction, respectively. u r and v r are the longitudinal and lateral speed through water, respectively. B ¼ arctan v=u ð Þis the sideslip angle. d denotes the rudder angle, where d j j 35 is the input saturation.
The ship's MMG model was proposed by Ship Manoeuvring Mathematical Model Group 8 in the late 1970s, which took the rudder angle as the control input. Considering the wind, ocean current and wave disturbances, as well as servo response characteristics, a ship's MMG model could be expressed as where m is the ship mass. m x , m y , I zz , and J zz are the corresponding added masses and added moment of inertia, respectively. X and Y denote the forces in x and y frames, N denotes the moment, and W denotes the wind. It is noted that only the slowly varying forces are counteracted by the steering system. The oscillatory motion caused by the first-order waveinduced forces should be prevented from entering the feedback loop. H, P, and R denote the low-frequency hydrodynamic forces on the hull, propeller, and rudder, respectively. d r is the command rudder angle, K E is the steering gear control gain, and T E is the servo time constant. And X R , Y R , and N R are the rudder forces (moments) given by where t R is the factor of rudder resistance deductions, a H is the ratio of the hull's additional lateral force to the rudder lateral force caused by steering, x H is the distance from the steering induced hull lateral force center to the ship's center of gravity, and F N is the rudder positive pressure. The winds (moments) X W , Y W , and N W in the dynamics (1) are given 29 by where r a is the air density, a R is the wind angle of attack, U R is the relative wind speed, and A f and A s are the frontal and lateral projected areas, respectively. L oa is the length overall, and C wx (a R ), C wy (a R ), and C wn (a R ) are the nondimensional wind coefficients, respectively.
The waves (moments) X Wave , Y Wave , and N Wave are given 30 by where l is the wavelength, is the encounter angle, r is the seawater density, a is the wave amplitude, L is the ship length, and C Xw (l), C Yw (l), and C Nw (l) are coefficients related to wave length l, respectively.

Control design and stability analysis
The structure of the control In this article, the path following controller is designed by MPC that consists of internal model, feedback correction, and dynamic optimization. Its main advantage is that the input constraint could be handled explicitly 26,31 because MPC considers a standard QP problem with saturation only. First of all, the future states over the prediction horizon N p can be predicted according to the internal model and current states of the system. Then, the cost function including prediction errors would be built via the future states. Finally, the optimal control input will be obtained by solving the cost function with constraint. At the next sampling time, the same procedure is repeated on the basis of the feedback of the new measurements from the system, so that is also called rolling optimization. Path following pertains to following a desired path independent of time, without imposing any restrictions on the temporal propagation along the path. Thus, the revolution of propeller has been set a constant in this article. Therefore, the control objective becomes that a suitable rudder angle d is designed to force the ship to follow the reference path at the certain propeller revolutions. Moreover, the MMG model is employed directly as the internal model for the prediction of states on the basis of Euler iterative algorithm. The structure of the control is shown in Figure 2.
As shown in Figure 2, the output of ship motion system is the current states at the current sampling time, and the unknowns including parameters and external disturbances are approximated by RBF NN on the basis of the previous states. Then the internal model would predict the future states of the system over the prediction horizon N p based on internal model (MMG) and current states. After that, the cost function will be established with respect to the prediction error e. Finally, the optimal rudder angle with saturation will be obtained by solving the cost function in the optimal solution module. At the next sampling time, the same procedure will be repeated upon the feedback of the new measurements from the ship motion system.
Assumptions. The assumptions in this article are as follows: I. The states of ship such as x, y, ', u, v, and r can be measured. II. The derivative with respect to time of historical states such as _ u,_ v, and _ r could be obtained. III. The unknowns d u , d v , and d r are bounded, namely, d u j j d u max , d v j j d v max , and jd r j d r max .
Remark. In assumption I, the states could be measured through the sensors or observers in the navigation system. 20 Then, assumption II can be satisfied by collecting the historical data. Assumption III is reasonable and appeared in previous study. 10

MPC controller
First of all, the prediction equations of the states may be employed according to the Euler iteration algorithm 32 and equation (1) xðk where T c is the sampling time.
and _ f d are the discrete values of equation (1), respectively  whered u ðkÞ,d v ðkÞ, andd r ðkÞ are the estimated values of the total unknowns including parameters uncertainties and environmental disturbances. At the current sampling time k, we can predict the future lateral displacements over the prediction horizon N p according to equation (5) as follows Then, the prediction errors e y k þ j ð Þ; j ¼ 1; 2; . . . ; N P could be calculated via equation (7) and the reference path y d At the sampling time k, the cost function can be established by the prediction errors where Q is the weighting matrix. The optimal control law d Ã k þ 1 ð Þwill be calculated by solving the QP problem and the partial differential equation @f Ã e y ; d À Á =@d ¼ 0 with the constraint d min d d max . At each time over the prediction horizon, in most of the previous MPC applications, all the corresponding future inputs d k þ1 ; d k þ2 ; . . . ; d k þN P are considered during predicting future states, whereas the first optimal input d Ã kþ1 only is chosen as the final control law to the system. However, we only consider a common future input over the prediction horizon, namely And the control law is adjusted online without the analytical expression because MPC depends on the modern computer. 26 The RBF NN is used to compensate the unknowns such asd u k ð Þ,d v k ð Þ, and d r k ð Þ that they are included in f Ã e y ; d À Á . N p , N c , and T c , namely prediction horizon, control horizon, and prediction sampling time, respectively, are control parameters to be tuned.

Radial basis function neural network
The gradient descent method of RBF NN is used to approximate the unknowns includingd u k ð Þ,d v k ð Þ, andd r k ð Þ on the basis of the previous information. RBF NN can approximate both the linear and nonlinear continuous functions over a compact set. 28 Takingd r k ð Þ as an example, the RBF NN method is given by where y n is output vector, and here isd r k ð Þ. n is the NN node number, n > 1, and w ¼ [w 1 , w 2 , . . . , w n ] is weight vector. The vector of RBF h is given by where x is the input vector, namely , c j is the center of the receptive field, and b j is the width of Gaussian function where a 2 0; 1 ð Þis momentum factor, the method to determine w, b, and c is given by where h 2 0; 1 ð Þ is learning rate and y h is previous information. w, b, and c could be adjusted online with gradient descent method. The acquisition processes ofd u k ð Þ,d v k ð Þ are the same as that ofd r k ð Þ.
Stability analysis Lemma 1. The continuous function will obtain minimum value if it has a point x where the derivative equals zero (g 0 (x) ¼ 0) and the left and right derivative symbols of this point are opposite (g 0 (x À D) < 0 and g 0 (x þ D) > 0). 33 From equations (1) to (6), we can predict the value of where F r is a positive constant without d, f r is the sum of other terms without d, and g is the rectification coefficient. g and b R 34 could be given by where C b is the square coefficient and B is the ship breadth. It is reasonable to suppose that the sideslip angle b is not very large and velocity is almost constant, so let b j j 25 , U 7:5 m=s, and r j j 0:02 rad=s here, then the derivative of equation (14) with respect to d yields From equations (5) and (6), the value of The derivative of equation (17) with respect to d, and substituting equation (16) into (17) yields The second equation of (5) can be reformulated as follows Þ . Combining equations (8) and (9), the cost function will be built as follows y d does not contain d, and the component of rudder force in u or v is smaller compared to angular rate r, 4 so the partial derivative of equation (20) with respect to d yields where z 1 , z 2 , and z 3 are as follows z 3 is always negative because d j j 35 . Then, substituting equations (14) and (17) into equation (22) yields where cos(d) is positive due to d j j 35 . The sign of q e is only determined by P, N p , or T c because the sideslip angle b is small and the heading is constant at current time. The sign of P will be plus as d ¼ À35 , and q e will increase with the increasement of N p or T c . If N p and T c are selected suitably, then z 1 > 0, z 2 > 0, and @f =@d < 0. The sign of P would be minus as d ¼ 35 , and q e will decrease with the increasement of N p or T c . If N p and T c are selected suitably, z 1 will be negative, whereas z 2 is positive, and @f =@d > 0. So there are suitable parameters N p and T c to satisfy equation (21) as follows The process of ship motion is continuous, so f is also a continuous function. According to equation (24) and the zero point theorem, there is an optimal rudder angle d* to satisfy @f/@d* ¼ 0. And it can fulfill @f/@d < 0 as d < d*, or @f/@d > 0 as d > d*. Thus, f will be the minimum at d* from lemma 1. Finally, f and e y will converge to zero gradually by obtaining the minimum of the f at each time during rolling optimization, then y(k) will also converge to y d (k) gradually.

Numerical simulation
To illustrate the effectiveness of the designed MPC scheme, the ship "Yulong" is taken for an example in the simulation. 8 The ship's dynamic is given by equation (1)   are 83 m, ' þ 135 À 30 sin(0.02t), and 3 m, respectively. The parameters of the controller are as follows: N p ¼ 18, N c ¼ 1, and T c ¼ 2.0 s are the parameters of MPC, initial c 0 ¼ À1 À 0:5 0 0:5 1; À 1 À 0:5 0 0: The simulation results are shown in Figures 3 to 6, where y, ', and d denote the path, heading, and rudder angle, respectively. The heading ' and rudder angle d have the time-varying deviation for the reason of the influence of the time-varying disturbances such as wave, wind, and ocean current. The rudder angle d is smooth with more stringent limitation on the rudder amplitude constraint within 10 . This illustrates that the ability of the MPC could handle the input saturation effectively. In Figure 3, the path y has the small overshoot less than the ship breadth. The maximum path tracking fluctuation is less than 5 m or less than 20% of the ship width on passage. In short, the control scheme could force the underactuated ship to follow the desired path accurately with uncertain parameters and environmental disturbances.

Conclusions
This article proposed an MPC design method, aiming at a path following control for underactuated surface vessels with input saturation, parameters uncertainties, and environmental disturbances. The choice of MPC for this application is primarily motivated by the need to explicitly handle input saturation. The Euler iterative algorithm reduced the amount of calculation of MPC. The total unknowns, including parameters uncertainties and environmental disturbances, were approximated by RBF NN to enhance the system robustness. The MMG model was employed as a predictive model to improve the control accuracy. The numerical simulation illustrated the effectiveness of designed controller.

Declaration of conflicting interests
The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.