Synchronization Control of Dynamic Positioning Ships Using Model Predictive Control

Abstract

Underway replenishment is essential for ships performing long-term missions at sea, which can be formulated into the problem of leader-tracking configuration. Not only the position and orientation but also the velocities are required to be controlled for ensuring the synchronization; additionally, the control inputs are constrained. On this basis, in this paper, a novel synchronization controller on account of model predictive control (MPC) for dynamic positioning (DP) ships is devised to achieve underway replenishment. Firstly, a novel synchronization controller based on MPC is devised to ensure the synchronization of not only the position and orientation but the velocities; furthermore, it is a beneficial solution for its advantages in handling the control input constraints ignored in most studies of underway replenishment. Secondly, a neurodynamic optimization system is applied to the implementation of MPC, which can improve the computational efficiency and shorten the simulation time. Thirdly, stability, frequently neglected by traditional MPC, is ensured by the means of adding a terminal cost function exported from the Lyapunov equation into the objective function. Finally, the effectiveness and advantages of the proposed control design are illustrated by extensive simulations.

1. Introduction

Dynamic positioning (DP) control has received intense research interest due to its wide applications in sea exploration such as offshore oil, underwater cable laying, and gas drilling. It is well known that DP control is autonomous, which can regulate the horizontal position and course of marine vessels either bound to a specific point, namely, absolute DP, or to a moving unit, that is relative DP [1,2]. Considerable solutions have been presented for the absolute DP in the past years [2,3,4,5,6,7,8,9,10]. Nevertheless, the relative DP, like underway replenishment, is still with challenges. Underway replenishment operation involves cargo transfer between two ships in transit, which is essential for ships performing long-term missions at sea to avoid returning to the base for replenishing storage or personnel [11,12,13]. It is normally formulated into the leader-tracking configuration [14,15,16,17,18,19,20]. Specifically, the leader vessel is controlled to sail on the steady velocity and course independently; however, the follower vessel requires to move upward alongside with a predefined offset to receive cargos. Apparently, for underway replenishment, not only the position and orientation are required to be controlled precisely, but the velocities must be synchronized simultaneously. There are six variables that needed to be controlled, namely longitudinal displacement, lateral displacement, yaw angle, surge velocity, sway velocity, and yaw rate, which lead to a multi-variable control problem. Additionally, it is worth noting that the constraint of the control input is exiting in the actual underway replenishment, but this problem is ignored in the above studies [11,12,13].

In previous years, model predictive control (MPC) has captured more attention owing to the formidable advantages of handling the multi-variable control problem with constraints [21]. The purpose of MPC is to forecast the future state in advance in virtue of an explicit model. Then, a control sequence can be calculated by means of handling an optimal problem on-line using current information of the system at every sampling moment. Control horizon, that is, the first of the above calculated resultant control sequence, is used for the plant, and the entire procedure is regarded as the next control intervals [21,22,23]. MPC is a promising solution for the underway replenishment of the DP ship for the following reasons: (1) the input constraints embedded in the actuators of the DP ship are directly considered in the design of control actions; (2) the multiple control variables are incorporated into one objective function where the optimization is conducted; (3) the negative effects of model uncertainty and disturbance can be compensated by the on-line calculation process; (4) the predictive capability of MPC is advantageous for sailing safety because the control action can be made in advance. The major challenges that prevent MPC from being widely applied in practice are computational efficiency and stability [21,22,23]. Fortunately, the neurodynamic optimization based on the parallel structure of a recurrent neural network proves to be efficient for the optimization problem [24,25,26], which has been applied to various areas successfully [26,27,28,29]; therefore, it is applied herein. For the sake of stability, in this paper, the MPC with a quasi-infinite horizon is introduced. Specifically, a linear state feedback gain is first calculated; then the terminal penalty matrix, which is the primary part of the terminal cost function, is calculated on account of the solution of a Lyapunov equation; finally, the resultant terminal cost function is applied into the objective function [30,31].

In summary, the primary contributions are summarized blew:

(1)

The underway replenishment of the DP ship is formulated into the leader-tracking configuration, and a novel model predictive control design (MPC) is presented; both the position and orientation are controlled to achieve the desired responses, and the velocities are controlled to be synchronized simultaneously.

(2)

The control input constraints, which are ignored in most studies of underway replenishment [11,12,13], are considered in the presented MPC design herein; the comparative research is carried out to validate the advantages.

(3)

The closed-loop stability, frequently neglected by traditional model predictive control, can be guaranteed herein by means of integrating the terminal cost function exported from the Lyapunov equation into the objective function.

This paper can be organized into the following sections. Concretely, Section 2 depicts preliminaries and problem formulation. Section 3 describes the design of the synchronization controller on the basis of MPC for DP ship in underway replenishment. Section 4 proposes the simulation results and detailed analysis. Section 5 summarizes the conclusions.

2. Preliminaries and Problem Formulation

2.1. The Mathematical Model of DP Ship

Normally, 3-DoF (degree-of-freedom) ship motion is considered for the maneuvering of marine surface vessels, namely surge, sway, and yaw. The kinematics of 3-DoF ship motion can be described as [32]

$$\dot{\eta }=J(\psi )\upsilon$$

(1)

here $\eta ={\left[x,y,\psi \right]}^T$ denotes the lengthwise shift, transverse shift, and course separately; $\upsilon ={[u,v,r]}^T$ represents the surge velocity, the sway velocity, and the yaw rate respectively; the rotation matrix is

$$J(\psi )=\left[\begin{array}{ccc}\cos \psi & -\sin \psi & 0\\ \sin \psi & \cos \psi & 0\\ 0& 0& 1\end{array}\right]$$

(2)

The kinetics of ship motion is normally defined as [33]

$$M\dot{\upsilon }=-C(\upsilon )\upsilon -D(\upsilon )\upsilon +{\tau }_T$$

(3)

where $M=\left[\begin{array}{ccc}m-X_{\dot{u}}& 0& 0\\ 0& m-Y_{\dot{v}}& m{x}_g-Y_{\dot{r}}\\ 0& m{x}_g-N_{\dot{v}}& I_z-N_{\dot{r}}\end{array}\right]$ is the inertia matrix; the damping matrix is $D(\upsilon )=\left[\begin{array}{ccc}-X_u-X_{\left|u\right|u}\left|u\right|-X_{uuu}{u}^2& 0& 0\\ 0& -Y_v-Y_{\left|v\right|v}\left|v\right|-Y_{\left|r\right|v}\left|r\right|& -Y_r-Y_{\left|v\right|r}\left|v\right|-Y_{\left|r\right|r}\left|r\right|\\ 0& -N_v-N_{\left|v\right|v}\left|v\right|-N_{\left|r\right|v}\left|r\right|& -N_r-N_{\left|v\right|r}\left|v\right|-N_{\left|r\right|r}\left|r\right|\end{array}\right]$, and $C(\upsilon )=\left[\begin{array}{ccc}0& 0& -(m-Y_{\dot{v}})v-(m{x}_g-Y_{\dot{r}})r\\ 0& 0& (m-X_{\dot{u}})u\\ (m-Y_{\dot{v}})v+(m{x}_g-Y_{\dot{r}})r& -(m-X_{\dot{u}})u& 0\end{array}\right]$ denotes the Coriolis and centrifugal matrix; furthermore, ${\tau }_T={[\begin{array}{ccc}{\tau }_1& {\tau }_2& {\tau }_3\end{array}]}^T$ are the control inputs, and it has $\left|{\tau }_1\right|\le 2.0\mathrm{N}$, $\left|{\tau }_2\right|\le 2.0\mathrm{N}$, and $\left|{\tau }_3\right|\le 1.5\mathrm{N}\cdot \mathrm{m}$; The $X_{(\cdot )}$, $Y_{(\cdot )}$, $N_{(\cdot )}$ represent hydrodynamic parameters [33].

2.2. Wave-Frequency (WF) Model

The WF model, which is widely adopted for stimulating the disturbances for dynamic positioning [2,32,34,35], is employed herein. The basic idea is to use a second-order system as an approximation for the linear wave spectrum for each degree of freedom. It has the following form [2,32]

$$\begin{array}{l}{\dot{\xi }}_x=x_w\\ {\dot{\xi }}_y=y_w\\ {\dot{\xi }}_{\psi }={\psi }_w\\ {\dot{x}}_w=-2\varsigma w_0{x}_w-w_0^2{\xi }_x+K_w{w}_{n1}\\ {\dot{y}}_w=-2\varsigma w_0{y}_w-w_0^2{\xi }_y+K_w{w}_{n2}\\ {\dot{\psi }}_w=-2\varsigma w_0{\psi }_w-w_0^2{\xi }_{\psi }+K_w{w}_{n3}\end{array}$$

(4)

where $w_0$ is the frequency of dominating wave, $\varsigma$ is the damping ratio of dynamic model, $w_{n1}$ is Gaussian noise, and $K_w$ is the parameter related to the wave intensity (the relationship between parameters’ value and sea state are given Table 4 in [2]). The disturbances can be presented as ${\eta }_w={\left[\begin{array}{ccc}{x}_w& y_w& {\psi }_w\end{array}\right]}^T$. Above all, the total motion of DP ship can be represented as [2,32]

$${\eta }_{total}=\eta +{\eta }_w$$

(5)

2.3. Problem Formulation

The underway replenishment control problem herein is formulated into a leader-tracking configuration [36,37]. The leader ship is controlled independently to maintain its course and speed, and the tracking ship is controlled to move up alongside a given offset. The purpose of the traditional tracking problem is to make the position of a marine ship follow a predefined trajectory as closely as possible; however, for missions, such as underway replenishment, the major challenge is to require the synchronization of velocities simultaneously while controlling the position and direction precisely. The leader ship is regarded as the reference for the follower ship, and the model can be depicted as

$$\begin{array}{c}{\dot{\eta }}_r=J({\psi }_r){\upsilon }_r\\ M{\dot{\upsilon }}_r=-C({\upsilon }_r){\upsilon }_r-D({\upsilon }_r){\upsilon }_r+{\tau }_{Tr}\end{array}$$

(6)

By following the coordinate transformation [16,38,39], Equation (1) can be rewritten as

$$\begin{array}{l}{z}_1=x\cos \psi +y\sin \psi \\ z_2=-x\sin \psi +y\cos \psi \\ z_3=\psi \end{array}$$

(7)

and the kinematics part of Equation (6) with the offset displacement $(x_o,y_o)$ can be represented as

$$\begin{array}{l}{z}_{1r}=(x_r-x_o)\cos {\psi }_r+(y_r-y_o)\sin {\psi }_r\\ z_{2r}=-(x_r-x_o)\sin {\psi }_r+(y_r-y_o)\cos {\psi }_r\\ z_{3r}={\psi }_r\end{array}$$

(8)

Then, the error system of synchronization control of the DP ship for underway replenishment can be defined as

$$\begin{array}{l}{z}_{ie}=z_i-z_{ir}i=1,2,3\\ u_e=u-u_r\\ v_e=v-v_r\\ r_e=r-r_r\end{array}$$

(9)

MPC is a natural solution for the multi-variable control problem of Equation (9) because the variables can be treated as an objective function systematically where the optimization can be applied [21,22,23,40,41,42]. Additionally, the control inputs of the DP ship are constrained in practice, so it is particularly important for sailing in proximity. Considering the large inertia of ship motion is the main cause for collision, the predictive mechanism of MPC allows it to handle the constraint ahead of time, which can mitigate the influence of the large inertia of ship motion. Therefore, in this paper, a novel synchronization controller on account of model predictive control (MPC) for dynamic positioning (DP) ships is devised to achieve underway replenishment. Firstly, a novel synchronization controller based on MPC is devised to ensure the synchronization of not only the position and orientation but the velocities; furthermore, it is a beneficial solution for its advantages in handling the control input constraints ignored in most studies of underway replenishment. Secondly, a neurodynamic optimization system is applied to the implementation of MPC, which can improve computational efficiency and shorten the simulation time. Thirdly, stability, frequently neglected by traditional MPC, is ensured by the means of adding a terminal cost function exported from the Lyapunov equation into an objective function. Additionally, the design block diagram of the presented synchronization controller is shown in Figure 1.

3. The Controller Design

In this part, the synchronization controller of the DP ship for underway replenishment is devised on the basis of the stable MPC and neurodynamic optimization.

By differentiating Equation (9) and linearizing the resultant equation, it yields

$${\dot{x}}_t=A_t{x}_t+B_t{\tau }_T$$

(10)

where $x_t={\left[\begin{array}{cccccc}{z}_{1e}& z_{2e}& z_{3e}& u_e& v_e& r_e\end{array}\right]}^T$, and $A_t=\left[\begin{array}{cccccc}0& 0& 0& 1& 0& 0\\ 0& 0& 0& 0& 1& 0\\ 0& 0& 0& 0& 0& 1\\ 0& 0& 0& -\frac{d_{11}}{m_{11}}& 0& 0\\ 0& 0& 0& 0& -\frac{d_{22}}{m_{22}}& 0\\ 0& 0& 0& 0& 0& -\frac{d_{33}}{m_{33}}\end{array}\right]$, $B_t=\left[\begin{array}{ccc}0& 0& 0\\ 0& 0& 0\\ 0& 0& 0\\ \frac{1}{m_{11}}& 0& 0\\ 0& \frac{1}{m_{22}}& 0\\ 0& 0& \frac{1}{m_{33}}\end{array}\right]$, here $m_{11}=m-X_{\dot{u}}$, $m_{22}=m-Y_{\dot{v}}$, $m_{33}=I_z-N_{\dot{r}}$, $d_{11}=-X_u$, $d_{22}=-Y_v$, $d_{33}=-N_r$. Based on a sampling moment, Equation (10) is rewritten as the following discrete form

$$x(k+1)=Ax(k)+B\tau (k)$$

(11)

where $x,y\in R^6$, $A\in R^{6\times 6}$, $B\in R^{6\times 3}$; $k$ represents the current sampling moment.

3.1. Design of Terminal State Matrix

A quadratic terminal cost function is devised to guarantee stability, and the terminal state matrix, which is the primary part of the terminal cost, is derived from an appropriate Lyapunov equation [30,31]. The process of design is described as follows.

A linear state feedback law for Equation (11) can be determined by

$$\tau =Kx$$

(12)

where $K$ is the vector of control gain. Then, the state matrix of terminal cost $P$ can be derived from the following Lyapunov equation

$${(A_k+\kappa I)}^{T}P+P(A_k+\kappa I)=-(Q_{lq}+K^T{R}_{lq}K)$$

(13)

where $A_k=A+BK$; $Q_{lq}$ and $R_{lq}$ are weighing matrices on control output and input [30,31]. $\kappa$ is chosen as according to the following principle

$$\kappa <-{\lambda }_{\max }(A_k)$$

(14)

It ensures the negativity of the real parts of all eigenvalues of $A_k+\kappa I$ [30,31].

3.2. MPC Design

Given that $N_P$ is the predictive horizon, $N_c$ is the control horizon [40,41,42]. Then, the further responses of Equation (11) can be forecasted by

$$\begin{array}{cc}\hfill x(k+1& |k)=Ax(k)+B\tau (k)\hfill \\ \hfill x(k+2& |k)=Ax(k+1|k)+B\tau (k+1)\hfill \\ & =A^{2}x(k)+AB\tau (k)+B\tau (k+1)\hfill \\ & ⋮\hfill \\ \hfill x(k+N_P& |k)=A^{N_P}x(k)+A^{N_P-1}B\tau (k)+A^{N_P-2}B\tau (k+1)\hfill \\ & +\dots +A^{N_P-N_c}B\tau (k+N_c-1)\hfill \end{array}$$

(15)

where $x(k+i|k),i=1,2,\cdots ,N_P$ is the forecast of $x(k+i)$.

The following vectors can be described as

$$\begin{array}{l}X={[\begin{array}{cccc}x{(k+1|k)}^T& x{(k+2|k)}^T& \cdots & x{(k+N_P|k)}^T\end{array}]}^T\in R^{6N_P}\\ \mathsf{\Gamma }={[\begin{array}{cccc}\tau {(k)}^T& \tau {(k+1)}^T& \cdots & \tau {(k+N_c-1)}^T\end{array}]}^T\in R^{3N_c}\end{array}$$

(16)

For the sake of simplicity, Equation (15) is rewritten as the following form

$$X=Fx(k)+\mathsf{\Phi }\mathsf{\Gamma }$$

(17)

where $F=\left[\begin{array}{c}A\\ A^2\\ ⋮\\ A^{N_p}\end{array}\right]\in R^{6N_P\times 6}$, $\mathsf{\Phi }=\left[\begin{array}{cccc}B& 0& \cdots & 0\\ AB& B& \cdots & 0\\ ⋮& \ddots & & ⋮\\ A^{N_p-1}B& A^{N_p-2}B& \dots & A^{N_p-N_c}B\end{array}\right]\in R^{6N_P\times 3N_c}$.

Then, the multi-variable control problem of the synchronization of the DP ship is transformed into the optimization problem with a terminal cost

$$J=X^T{Q}_{mpc}X+{\mathsf{\Gamma }}^T{R}_{mpc}\mathsf{\Gamma }+X_P^{T}P{X}_P$$

(18)

subject to

$${\mathsf{\Gamma }}_{\min }\le \mathsf{\Gamma }\le {\mathsf{\Gamma }}_{\max }$$

(19)

where $Q_{mpc}$, $R_{mpc}$ are weighing matrices, $X_P^{T}P{X}_P$ denotes the terminal cost, $X_P$ is the terminal state, $P$ is the state matrix of terminal cost determined by Equation (13).

3.3. Neurodynamic Optimization Design

One challenge of the implementation of the MPC method is the computational burden. In this section, a neurodynamic system is set up to deal with the constrained optimization problem parallelly using the PNN.

Substituting Equation (17) into Equation (18), it yields

$$\begin{array}{cc}\hfill J& ={\left[Fx(k)+\mathsf{\Phi }\mathsf{\Gamma }\right]}^T{Q}_{mpc}\left[Fx(k)+\mathsf{\Phi }\mathsf{\Gamma }\right]+{\mathsf{\Gamma }}^T{R}_{mpc}\mathsf{\Gamma }\hfill \\ & +{\left[F{x}_p(k)+\mathsf{\Phi }\mathsf{\Gamma }\right]}^{T}P\left[F{x}_p(k)+\mathsf{\Phi }\mathsf{\Gamma }\right]\hfill \\ & ={\left[Fx(k)\right]}^T{Q}_{mpc}\left[Fx(k)\right]+{\mathsf{\Gamma }}^T{\mathsf{\Phi }}^T{Q}_{mpc}\mathsf{\Phi }\mathsf{\Gamma }+{\mathsf{\Gamma }}^T{R}_{mpc}\mathsf{\Gamma }+{\left[F{x}_p(k)\right]}^{T}P\left[F{x}_p(k)\right]\hfill \\ & +{\mathsf{\Gamma }}^T{\mathsf{\Phi }}^{T}P\mathsf{\Phi }\mathsf{\Gamma }+2{\left[Fx(k)\right]}^T{Q}_{mpc}\mathsf{\Phi }\mathsf{\Gamma }+2{\left[F{x}_p(k)\right]}^{T}P\mathsf{\Phi }\mathsf{\Gamma }\hfill \\ & ={\mathsf{\Gamma }}^T({\mathsf{\Phi }}^{T}P\mathsf{\Phi }+{\mathsf{\Phi }}^T{Q}_{mpc}\mathsf{\Phi }+R_{mpc})\mathsf{\Gamma }+\left[2{\left[F{x}_p(k)\right]}^{T}P\mathsf{\Phi }+2{\left[Fx(k)\right]}^T{Q}_{mpc}\mathsf{\Phi }\right]\mathsf{\Gamma }\hfill \\ & +{\left[Fx(k)\right]}^T{Q}_{mpc}\left[Fx(k)\right]+{\left[F{x}_p(k)\right]}^{T}P\left[F{x}_p(k)\right]\hfill \end{array}$$

(20)

Define

$$\begin{array}{l}H={\mathsf{\Phi }}^{T}P\mathsf{\Phi }+{\mathsf{\Phi }}^T{Q}_{mpc}\mathsf{\Phi }+R_{mpc}\\ W=2{\left[F{x}_p(k)\right]}^{T}P\mathsf{\Phi }+2{\left[Fx(k)\right]}^T{Q}_{mpc}\mathsf{\Phi }\\ \xi ={\left[Fx(k)\right]}^T{Q}_{mpc}\left[Fx(k)\right]+{\left[F{x}_p(k)\right]}^{T}P\left[F{x}_p(k)\right]\end{array}$$

(21)

Thus, Equation (20) can be written as

$$J={\mathsf{\Gamma }}^{T}H\mathsf{\Gamma }+W\mathsf{\Gamma }+\xi$$

(22)

The control variable $\mathsf{\Gamma }$ is considered as the states of the PNN, and the neurodynamic system is described as

$$\frac{d\mathsf{\Gamma }}{dt}=\vartheta \left\{P_{\mathsf{\Omega }}(\mathsf{\Gamma }-\mathsf{\Omega }(\mathsf{\Gamma }))-\mathsf{\Gamma }\right\}$$

(23)

here $\vartheta$ is the convergent rate, $\mathsf{\Omega }(\mathsf{\Gamma })$ is set as $\nabla J$ [24,25], and $P_{\mathsf{\Omega }}$ is the projection operator.

$$P_{\mathsf{\Omega }}({\tau }_i)=\{\begin{array}{l}{\tau }_{\min },{\tau }_i<{\tau }_{\min }\\ {\tau }_i,{\tau }_{\min }\le {\tau }_i\le {\tau }_{\max }\\ {\tau }_{\max },{\tau }_i>{\tau }_{\max }\end{array}$$

(24)

Consider the following Lyapunov function

$$V(\mathsf{\Gamma })={\displaystyle {\int }_0^1(\mathsf{\Gamma }-{\mathsf{\Gamma }}^{*})}\mathsf{\Omega }({\mathsf{\Gamma }}^{*}+t(\mathsf{\Gamma }-{\mathsf{\Gamma }}^{*}))dt$$

(25)

By differentiating Equation (25), we can obtain [24,25],

$$\frac{dV(\mathsf{\Gamma })}{d\mathsf{\Gamma }}=\nabla V{(\mathsf{\Gamma })}^T=\mathsf{\Omega }$$

(26)

where $\nabla V(\mathsf{\Gamma })$ represents the gradient of $V$.

Then, we have

$$\frac{dV(\mathsf{\Gamma })}{dt}=\left[\frac{dV}{d\mathsf{\Gamma }}\right]\frac{d\mathsf{\Gamma }}{dt}=\vartheta \mathsf{\Omega }{(\mathsf{\Gamma })}^T\left\{P_{\mathsf{\Omega }}(\mathsf{\Gamma }-\mathsf{\Omega }(\mathsf{\Gamma }))-\mathsf{\Gamma }\right\}$$

(27)

Since $\mathsf{\Omega }{(\mathsf{\Gamma })}^T\left\{P_{\mathsf{\Omega }}(\mathsf{\Gamma }-\mathsf{\Omega }(\mathsf{\Gamma }))-\mathsf{\Gamma }\right\}\le -\vartheta {\Vert P_{\mathsf{\Omega }}(\mathsf{\Gamma }-\mathsf{\Omega }(\mathsf{\Gamma }))-\mathsf{\Gamma }\Vert }^2$ [24,25]. Thus,

$$\frac{dV}{dt}\le 0$$

(28)

We can conclude that the neurodynamic system is steady. In other words, each trajectory corresponding to the dynamic system can converge to an equilibrium point from any initial point in the domain [24,25,26,27,28,29].

4. Simulation Results

Multiple simulations are carried out to illustrate the effectiveness and advantages of the presented novel synchronization controller on account of model predictive control. Specifically, the multiple simulations can be segmented into two groups: the first group of simulations is to illustrate the effectiveness of the presented synchronization controller from the different initial points with disturbances, the disparate predictive horizon, and the changes in velocities and course; the advantages of the presented synchronization con-troller are illustrated by the two comparative simulations, the first is between PID, LQR, and the presented synchronization controller with and without input constraints, and the other is between conventional MPC and the presented MPC integrating with PNN. The parameters of the simulation model are given in [33], specifically, $T_s=0.01s$, $N_c=5$, $N_P=180$, $Q_{mpc}=diag\left\{\begin{array}{cccccc}3000& 3000& 3000& 1& 1& 1\end{array}\right\}$, and $R_{mpc}=diag\left\{\begin{array}{ccc}1& 1& 1\end{array}\right\}$. And the parameter of PNN is set as $\vartheta =0.001$.

4.1. Performance of Synchronization Control from Different Initial Points with Disturbances

In this section, the control performance of the presented novel synchronization controller based on MPC starting from different initial points with disturbances is studied to illustrate the effectiveness of the presented synchronization controller. The initial states of leader and follower are given in

Table 1. Initial states of leader and follower.

	Leader Ship	Follower Point-1	Follower Point-2	Follower Point-3
$x$ (m)	0	18	25	30
$y$ (m)	0	0	−5	−10
$\psi$ (deg)	45	45	45	45
$u$ (m/s)	0.1	0.5	1	1
$v$ (m/s)	0.1	0.1	0	1
$r$ (rad/s)	0.05	0.1	0	0.5

.

It is worth noting that the disturbances are simulated using the WF model in Section 2. The Beaufort number (BN) is chosen as 5, and the environmental effects of BN 5 can be described with wind speed 10.7 m/s, significant wave height $H_s=2.1\mathrm{m}$, peak wave period $T_0=7.5\mathrm{s}$, and current speed 0.75 m/s [43]. The corresponding parameters of WF model are therefore set as ${\omega }_0=0.838$, $\varsigma =0.257$, and $K_w=0.296$ (the values of parameters for different Beaufort numbers are given in [2]).

In addition, the simulation results are shown in Figure 2, Figure 3 and Figure 4. Specifically, Figure 2 denotes the tracking performance of the presented synchronization control from different initial points with disturbances. Figure 3 represents the synchronization of velocities. Additionally, Figure 4 indicates the corresponding responses of control inputs. For traditional underway replenishment, normally, it involves two separated stages. First, the follower ship is controlled to the predefined location in the formation independently. Then, the velocities are controlled to synchronize with the leader ship. In this section, the position and velocities of the follower ship are designed to be control outputs together. As shown in Figure 2 and Figure 3, the follower vessel can move along with the leader ship with a given offset from different initial points with disturbances automatically. Meanwhile, the synchronization of velocities is achieved. After all, the effectiveness of the presented synchronization controller is illustrated by the simulation results.

4.2. Influence of Prediction Horizon

The predictive horizon $N_P$ is the key control parameter of MPC, representing the length of prediction. The ship motion normally suffers from large inertia, which is rather dangerous for navigation safety. The predictive capability of MPC makes it a favorable solution, because the control action can be decided ahead of time by virtue of the predicted further responses. This is crucial for underway replenishment to avoid the collision. Therefore, in this section, the simulations are carried out by set different predictive horizon $N_P$ varies as 100, 400, and 450 to illustrate the effectiveness of the presented synchronization controller. It is worth noting that $N_P=100$ indicates that the future responses in 1 s ($T_s=0.01\mathrm{s}$) are predicted and then used in the control design. The simulation results are shown in Figure 5, Figure 6 and Figure 7, which are the tracking performance, synchronization of velocities, and the corresponding responses of control inputs with different predictive horizons.

It is worth noting that the control input constraints $\left|{\tau }_1\right|\le 2.0\mathrm{N}$, $\left|{\tau }_2\right|\le 2.0\mathrm{N}$, and $\left|{\tau }_3\right|\le 1.5\mathrm{N}\cdot \mathrm{m}$ are considered in this case and the forthcoming cases. As we can see, when $N_P=100$, the responses of the follower ship fluctuate, particularly in the surge velocity, see Figure 6; when $N_P=400$ and $N_P=450$, the responses of the follower ship are nearly the same. It can be concluded that when $N_P\ge 400$, the control performance does not improve much as $N_P$ increases. Therefore, in this paper, the predictive horizon is set as $N_P=400$, i.e., the control action can be taken in advance approximately 4 s.

Note that, in Section 4.1, the predictive horizon is $N_P=180$ and the control performance seems fine because the control input constraints are not initiated, see Figure 4. The authors use such arrangement to show the influences of control input constraints, which is frequently ignored in underway replenishment [11,12,13]. However, it truly exists in reality. The details are reported the following case.

4.3. Comparative Studies on Input Constraints

In this section, the comparative simulations are carried out between proportion-integration-differentiation (PID) [11,12], linear quadratic regulation (LQR), and the presented synchronization controller to illustrate the advantages of the presented synchronization controller in solving input constraints. The parameters of PID are set as $K_{P1}=1000$, $K_{i1}=0.1$, $K_{d1}=10$, $K_{P2}=200$, $K_{i2}=0.1$, $K_{d2}=10$, $K_{P3}=200$, $K_{i3}=0.1$, and $K_{d3}=10$. The parameters of LQR are $Q=diag\left\{\begin{array}{cccccc}1000& 1000& 1000& 1& 1& 1\end{array}\right\}$, $R=diag\left\{\begin{array}{ccc}1& 1& 1\end{array}\right\}$.

Two scenarios, i.e., scenario (a): comparative studies without input constraints and scenario (b): comparative studies with input constraints, are conducted. The comparative performance is shown in Figure 8, Figure 9 and Figure 10. In scenario (a), the control performance is just slightly different using PID, LQR, and the presented control design. However, when the input constraints are initiated, the performance are very different, particularly in the synchronization of velocities, see Figure 9. Note that the length of simulations is 100 s in this section. However, in order to show the results clearly, only the transit process, where the main difference appears, is given in Figure 9 and Figure 10.

In Figure 9 scenario (b), as we can see, from the perspective of settling time and overshoot, the best control performance is realized by means of the presented MPC design. The reason is that the concept of MPC is to use the predicted future responses to determine current control action, and in this way, the violation of constraints can be predicted and avoided in advance. The LQR is better than PID because it is an optimal control method that can handle the multi-variable control problem such as underway replenishment in a way that PID does not attempt. For underway replenishment, not only the position and orientation but also velocities are required to be synchronized; moreover, the control input constraints exist in practice. Nevertheless, the existing references [11,12] use the PID method to design the controller, while control input constraints are not considered in all references [11,12,13]. It motivates the present study.

4.4. Performance of Synchronization Control with Changes in Velocities and Course

In this section, the control performance of the presented novel synchronization controller based on MPC with changes in velocities and course is studied to illustrate the effectiveness of the presented synchronization controller. The reference state information from the leader vessel is set up as follows

$$\{\begin{array}{cc}{u}_r=1\mathrm{m}/\mathrm{s},v_r=0\mathrm{m}/\mathrm{s},r_r=0\mathrm{r}\mathrm{a}\mathrm{d}/\mathrm{s}& 1\le t<50\mathrm{s}\\ u_r=0.5\mathrm{m}/\mathrm{s},v_r=0\mathrm{m}/\mathrm{s},r_r=-0.005\mathrm{r}\mathrm{a}\mathrm{d}/\mathrm{s}& 50\le t<100\mathrm{s}\\ u_r=1\mathrm{m}/\mathrm{s},v_r=0\mathrm{m}/\mathrm{s},r_r=0\mathrm{r}\mathrm{a}\mathrm{d}/\mathrm{s}& 100\le t<150\mathrm{s}\end{array}$$

(29)

The simulation results are described in Figure 11, Figure 12 and Figure 13. Specifically, Figure 11 represents the tracking performance of synchronization control with changes in velocities and course. Figure 12 shows the synchronization of velocities with changes in velocity and course. Figure 13 shows the corresponding responses of the control inputs. As we can see that the purpose of synchronization can be realized even the course and velocities change over time, which can illustrate the effectiveness of the presented synchronization controller.

4.5. Comparative Studies on Computation Efficiency

The computation efficiency is important for real-world applications. Therefore, in this section, the contrastive research between conventional MPC and the presented synchronization controller is carried out to illustrate the advantages of the presented synchronization controller in improving the computation efficiency (Both cases with Intel Core (TM) i5-6300HQ CPU, 2.3 GHz, 4 GB). It should be noted that the predictive horizon $N_P=400$ is the same as the above cases, as well as $N_c=400$.

Figure 14 and Figure 15 show the contrastive performance and control inputs of conventional MPC and the presented MPC integrating with PNN methods. As we can see, the simulation time of the presented MPC with the PNN method is much less than the traditional MPC, whereas the control performance of the two methods is nearly the same. It can be concluded that the presented MPC with the PNN method is with advantage in computational efficiency.

5. Conclusions

A novel model predictive control for underway replenishment of the DP ship has been presented in this paper. The underway replenishment has been formulated into a leader–follower configuration. A quasi-infinite horizon technique has been employed to guarantee stability by designing an appropriate terminal cost function based on the Lyapunov theorem. The simulation results have shown that the follower can move along with the leader ship automatically from different initial points under disturbances, and the synchronization of velocities can be achieved; then, the key parameter of the MPC has been investigated. Additionally, the contrastive research has illustrated the superiorities of the presented method in dealing with input constraints and computational intensity. Nevertheless, the presented method is limited to underway replenishment between the homogeneous ship model. In future, the further research will concentrate on applying the presented synchronization control method to achieve the underway replenishment between the heterogeneous ship model.