A Receding Horizon Navigation and Control System for Autonomous Merchant Ships: Reducing Fuel Costs and Carbon Emissions under the Premise of Safety

Abstract

In order to solve the multi-objective planning and trajectory tracking control problem related to maritime autonomous surface ships (MASSs), a new design scheme for autonomous navigation is proposed in this paper, with a receding horizon navigation and control (RHNC) system that contains navigation and control modules. In the navigation module, we designed a superposition field gradient descent local search algorithm based on the cost field, emission field, and guidance field to navigate the MASS reference path, and in the control module, we designed a nonlinear controller that can handle multiple constraints based on the NMPC framework. Under the new scheme, the navigation module completes local path planning to reduce costs and emissions, the control module accomplishes accurate trajectory tracking and real-time collision avoidance, and the information is transmitted in both directions between the two modules to collaboratively complete the MASS navigation and control tasks. We conducted a simulation study of the navigation algorithm and controller and the autonomous navigation system using a Kriso Container ship (KCS). The simulation results demonstrate the effectiveness of the proposed cooperative design scheme in reducing navigation costs and emissions and avoiding autonomous collision avoidances.

1. Introduction

Reducing costs, emissions, and collision avoidance in ship navigation has always been an hot topic in the field of marine engineering research, and the regulations on energy consumption and emissions issued by the International Maritime Organization (IMO), such as the Energy Efficiency Existing ship Index (EEXI) and the Carbon Intensity Indicator (CII), have brought the costs and emissions incurred in ship navigation into focus. With the development of artificial intelligence, unmanned vehicles, and advanced communication technologies, these hot issues are also at the forefront of research on maritime autonomous surface ships (MASSs). To achieve the goal of economic, green, and safe navigation, the route design becomes a complex multi-objective task, which is very challenging for the MASS navigation system. In order to reduce this burden, more research is needed, for example, by putting the safety objectives of navigation into the control system to solve. Therefore, the navigation and control modules in the MASS navigation system need to be designed to accomplish multi-objective autonomous navigation.

Guidance, navigation, and control (GNC) [1], the classical structure in autonomous navigation systems, contains a navigation module and a control module. Under the GNC structure, the path planning algorithm of the navigation module determines the reference path for the generation of trajectory information, which then determines the design of the controller and controls the control decision variables by the designed controller to finally maneuver the ship and track the reference trajectory. Path planning and trajectory tracking are the underlying processes of the navigation and control modules, which at this stage have generated abundant research on path planning and tracking control.

Path planning is usually divided into global and local path planning based on the type of environmental information [2]. Global path planning finds a safe path from the initial state to the target state considering known obstacles and assuming a complete environment model is available [3], so it is widely used in guidance systems for unmanned ships, robots, autonomous underwater vehicles, and unmanned vehicles. Grid-based global path planning was originally proposed in the form of Dijkstra’s algorithm [4,5]. The latter is a heuristic algorithm compared to Dijkstra’s algorithm, which also incorporates obstacle avoidance. Refs. [6,7] used the A* algorithm to plan weather routes, and in [8], it was applied to the path planning of unmanned surface vehicles (USVs). A constrained A* method was proposed for the optimal path planning of unmanned vehicles in constrained marine environments, with a vector map representing sea currents superimposed on a binary map to simulate static environmental disturbances, which are robust and computationally efficient for the global path planning of unmanned vessels under moving obstacles. An APF-based multi-directional a* algorithm was proposed in [9] to generate safer and shorter navigation trajectories within wind farm waters. To address the limitations of existing global path planning methods in which path planning is independent of the USV control phase, ref. [10] proposed an A*-VVGO algorithm that implements velocity variation by adding a temporal dimension to the map modeling process, which provides a new way to avoid obstacles compared to traditional path planning methods. In addition, in order to ensure the security of the path, the above algorithm introduces an artificial potential field method in map modeling.

From the above, it is clear that the A* algorithm has many applications and has been the subject of in-depth studies in the field of global path planning. In addition to the above two methods, the ant colony algorithm, which is an intelligent swarm algorithm, is also used for global path planning, which has good convergence and can find the globally optimal path [2,11,12]. The global path planning method is to find a suitable path between start and end points according to a search algorithm under the premise that the environmental information is known, which has limitations in scenarios where the environmental information is partially unknown and often requires improvements to achieve dynamic planning. In terms of safety, a static and dynamic obstacle layout in the environment is needed for global path planning before effective collision avoidance can be achieved.

Local planning methods often show good applicability when a ship only has access to local environmental information during navigation and needs to focus on the surrounding environment to avoid collisions with other ships. Local path planning algorithms include dynamic window [13], velocity obstacle [14], and artificial potential field (APF) [15]. The dynamic window algorithm (DWA) was used in [3] for local path planning to track the global path by tracking the intersection of global and local paths, avoiding dynamic obstacles. Artificial potential field (APF) methods are mathematically elegant and simple for dynamically complex conditions, but easily fall into local optimality; in [16], APF was improved so that it could be moved out of local optimality. A hybrid method was proposed in [17] based on the global planner A* algorithm and local planner APF method, called the multiple subgoal artificial potential field (MTAPF) algorithm, which can guarantee the optimality, rationality, and path continuity of the formation trajectory. This combined global and local planning method had significant optimization effects. A field simulating electromagnetic wave propagation was constructed in [18] via the fast marching method, which had no local minima problem and reduced the computational burden. Local path planning utilizes local environmental information around a ship, considers sensor information for situational awareness, and emphasizes avoiding nearby dynamic obstacles to generate feasible and safe paths. Often, the purpose of local path planning in existing studies is to avoid obstacles, especially dynamic obstacles.

Including carbon emissions under the premise of considering cost will inevitably make the ship’s routes more complicated. The challenges brought by carbon emissions coupled with ship control multivariate characteristics require us to use advanced nonlinear multi-constraint control methods as technical support for the control module to accomplish accurate trajectory tracking. In tracking control, most traditional research has focused on linear control techniques, such as PID [1]; however, ship dynamics are highly nonlinear, and nonlinear factors, such as wind and wave disturbances, and hydrodynamic damping can adversely affect the performance of the ship. At present, the research on ship motion control is mainly focused on the underdriven ship trajectory tracking control problem, which is better applied in the design of automatic rudders and track rudders. The commonly used robust control methods include sliding mode control (SMC), H∞ control, neural network control, and optimal control, among others [19]. A controller combining SMC and backstepping techniques was used in [20,21], showing good robustness. Ref. [22] designed a data-driven neural predictor based on integral parallel learning that does not require the support of a kinematic model. Ref. [23] proposed a method for robust neural network control with a new performance function, two additional terms, the radial basis function (RBF) neural network (NN), and the command filter backpropagation method, which can handle AUV trajectory tracking in unknown environments compared to conventional methods.

Although all of the above methods have good robustness and some of the neural network control algorithms can get rid of the dependence of traditional control algorithms on environmental information, a ship is subject to dynamic state constraints and the control constraints of rudder angle and main engine speed during motion, and there are certain limitations in the processing of these constraints via these methods. In contrast, model predictive control (MPC) can handle complex constraints very well [24,25]. In [26], the performance of H∞ and MPC was compared, and the conclusion was that the MPC controller performed better in terms of accuracy and response time under good driving conditions. For the nonlinear problem (NLP) of underactuated surface vessel tracking control, combined nonlinear model predictive control (NMPC) was used in [27] to transform static and dynamic object collision avoidance into a constraint problem embedded in MPC, further demonstrating the advantage of MPC in handling multiple constraints. An improved MPC was proposed in [28] based on global heading constraints (GCCs) and the event triggering mechanism (ETM), which ultimately optimized the deviation according to its predictive capability and energy saving-effect. Ref. [29] proposed a kind of collision avoidance constraint based on the ship safety domain, and a collision risk index (CRI) event trigger was added based on DCPA/TCPA, which is in turn based on the embedded collision avoidance constraint. This not only enables collision avoidance for dynamic and static obstacles of a smart ship but also reduces the computational complexity of the NMPC controller and the computation time of the system. The above studies illustrate that MPC, a control method that is more capable of handling multiple dynamics constraints, is well suited for MASS motion control.

In a study on combining navigation and control, the artificial potential field method and multi-constrained model predictive control (MMPC) were combined in [30] to complete path planning and trajectory tracking control, and it describes the tracking task as a multi-constraint model predictive control (MMPC) problem. In [31], the A* algorithm was combined with the sliding mode control approach as the underlying technique for an integrated path planning and trajectory tracking control framework. The path planning algorithm was used to find an optimal and smooth path while avoiding static and dynamic obstacles, and a pure path tracking approach was used to convert the generated geometric path into a time-dependent trajectory profile, which can serve as instructions for the trajectory tracking controller. In [32], the path planning problem was transformed into a rolling optimization problem. Combined with the NMPC controller, a unified rolling time-domain optimization scheme was proposed to solve the underwater survey path planning and trajectory tracking problem for AUVs. It is thus clear that the combination of navigation and control is ultimately a collaborative optimization of path planning and tracking control algorithms.

Most of the existing studies on path planning focus on global path planning with the goal of optimizing time [33] or energy consumption [34], or on safety-oriented local planning [35]. Studies also use a hybrid approach considering both [17], where global planning is mostly static, and local planning situations are mainly dynamic and static with collision avoidance planning, where safety is the top priority when avoiding collisions. As for the problem studied in this paper, the global and local planning algorithms have certain limitations. It is difficult to balance the three goals of cost reduction, emission reduction, and collision avoidance, and the combination of traditional control and path planning algorithms is not sufficient enough; there is only one-way information transfer between the two modules, which is not conducive to autonomous MASS navigation. Therefore, in this paper, the multi-objective path planning problem in MASS navigation is transformed into a search problem of the lowest point of the superposition potential field, the tracking control problem is transformed into an optimization problem of error dynamics in the design of the controller, and the ship dynamic and collision avoidance constraints are systematically added to this optimization problem, in order to take into account the three major objectives of ship navigation. The output of tracking control optimization is fed back into the path searching problem as information to provide starting state updates as the next moment of the path seeking optimization, forming a closed-loop structure between path planning and controller for collaborative optimization. In addition, in this paper, we also design an event trigger mechanism to close the potential cost and carbon emission field when the obstacle reaches the emergency hazard collision avoidance range of the ship, with navigation safety as the top priority.

The contributions of this study are as follows:

• The design of a new autonomous navigation system is proposed, which accomplishes the economic and green goals of navigation in the path planning module and ensures the safety of navigation in the control module;
• Based on the superposition of fields, cost and emission fields according to their consumption rate are constructed, and guided fields are designed to transform the path search problem into a gradient descent optimization problem to solve the multi-objective path planning problem of MASS;
• A navigation and control system framework for path planning and rolling NMPC tracking control co-optimization is designed for the problem of determining whether the reference trajectory and tracking controller track with high accuracy during multi-objective planning.

The rest of the paper is organized as follows: The problem statement is given in Section 2. Section 3 gives the specific procedure of establishing the potential field and the design of the path planning algorithm. Section 4 gives the design process of the three-degree-of-freedom model and controller for ship motion and collision avoidance constraints. Section 5 presents the design process of the real-time navigation and control system. Section 6 verifies the local path planning algorithm and the effectiveness of the system through simulation. Finally, Section 7 gives the conclusion.

2. Problem Statement

2.1. Problem Description

In the ocean, the density of ships is low, weather conditions are more complex compared with nearshore disturbance, and collision avoidance behavior is less frequent as well, so when a MASS is sailing in the ocean, most of the time, it is important to consider reducing fuel consumption and emissions during the voyage. Therefore, in this paper, we proposed a new solution to the MASS navigation and control problem, in which we transformed the multi-objective optimization problem into a superposition field seeking problem and put the collision avoidance problem into the control module to solve, which is closer to the actual movement of the ship than the navigation module. In addition, the ship’s motion system has the characteristics of time delay, large inertia, nonlinearity, and underdrive; therefore, precise tracking control of the complex trajectory of the ship has high requirements for both controller and path planning. It is often difficult for offline path planning to meet the precise trajectory tracking, and the synergy between navigation and control must be exploited to better accomplish the autonomous navigation task. Therefore, how to design a navigation and control system so that the ship can accurately track complex routes under multi-objective navigation decisions is the main problem to be solved in this paper.

From the above, it is clear that the mission of the MASS navigation system involves the following three aspects:

(i)

A guidance module for multi-objective path planning;

(ii)

A control module for precise trajectory tracking and to guarantee navigation safety;

(iii)

The collaborative optimization of guidance and control for information interaction and synergy optimization between navigation and control modules.

2.2. Framework of Receding Horizon Navigation and Control System

The navigation and control system framework proposed in this paper aimed toward multi-objective planning and trajectory tracking control for MASS ocean voyages and keeping the ship as close to the planned trajectory as possible while ensuring safety. At the beginning of each sampling moment, new planning was performed according to the change in the positional state of the ship to achieve real-time planning. Figure 1 gives the general architecture, and the core of the system is the navigation and the control modules.

In the navigation module, we transformed the wind, wave, and current environment of navigation into a superposed potential field covering costs and emissions as well as the shortest circuit target, performed a path search based on that field, turned the path points into trajectory points by adding time constraints, and transformed the reference trajectory points into state volume references via Equation (19). In the control module, we designed an objective optimization function for MASS trajectory tracking which contained two nonlinear performance metrics: the error between the reference state volume and the ship’s predicted state volume, and the rate of change in the control variables (rudder angle and propeller speed). The multiple constraints for nonlinear optimization included control variable constraints (rudder angle and propeller speed) and ship dynamic constraints [29]. In addition, to ensure that the path planning of multiple targets did not negatively affect autonomous collision avoidance, the collision avoidance mechanism was triggered when the radar detected an obstacle nearby, acting simultaneously on both the navigation and control modules: the cost potential field was turned off, and the collision avoidance constraint was turned on.

3. Guidance Module

Field theory is a concept of physics that involves the distribution of and changes in physical quantities in space [36], and the ability to superpose fields is applicable to the model construction of multiple targets; therefore, we designed the cost, emission, and guidance fields in this section, where the cost field and emission field are functions of the economic cost and increased CO${}_2$ emissions per unit time, and the guidance field is a function related to the distance of the target point and the “yaw” direction. The guidance field is a function related to the distance of the target point and the yaw direction (relative to the ship heading along the shortest path). The path search of the navigation module is based on the superimposed field environment, and the process of multi-objective route planning is the gradient descent process of the superimposed field. In addition, in order to prevent the path navigated by this module from affecting the collision avoidance performance of the subsequent control module when considering cost and emission reduction, we designed an event-triggering mechanism to switch between the regular navigation state and the collision avoidance state, and only the guidance field was retained to function in the latter.

3.1. Fields Based on Cost and Emission

The cost and emission fields are based on equations by [37,38] for the energy consumption of a ship while sailing at sea as a function of sea state, as follows (descriptions of parameters are given in

Table 1. Descriptions of parameters in equations.

Parameters	Meanings
K	Propeller efficiency
$\eta$	Propulsion efficiency (0.6 to 0.7)
$H_{\frac{1}{3}}$	Significant wave height
$P_s$	Still water power
$\rho$	Water density
$C_{ts}$	Still water drag coefficient
S	Wetted surface
$P_w$	Additional power for wave resistance
$C_w$	Drag coefficient for wave resistance
g	Vertical force
B	Ship width
L	Ship length
$P_a$	Additional power for wave resistance
$C_a$	Drag coefficient for aerodynamic
A	Surface area projected for wind
u	Wave speed
$u_a$	Wind speed
m	Cargo on board the vessel
n	Cargo weight constant
r	Wave constant (0.1)

):

$$\begin{array}{cc}\hfill & P=K(P_s+P_w+P_a)+P_{aux}\hfill \\ \hfill & K=\eta (v,H_{\frac{1}{3}})=max(\frac{1}{\eta (j+k·\sqrt{\frac{v}{v_d}})},\frac{1}{\eta (1-r·H_{\frac{1}{3}})})\hfill \\ \hfill & P_s=\frac{\rho ·C_{ts}·S·v^3}{2}(\frac{m·n}{DWT}+(1-n))\hfill \\ \hfill & P_w=\frac{C_w·\rho ·g·{(\frac{H_{\frac{1}{3}}}{2})}^2·B^2}{L}·(v+u)\hfill \\ \hfill & P_a=\frac{C_a·{\rho }_a·A·{(v+u_a)}^3}{2}\hfill \\ \hfill \end{array}$$

(1)

With known energy consumption (EC), the cost per unit time and carbon emissions (EM) can be determined by the following equation:

$$\begin{array}{cc}\hfill & EC=P·k_f·C_{fuel}+\frac{C_c}{24v_s}\hfill \\ \hfill \end{array}$$

(2)

$$\begin{array}{cc}\hfill & EM=P·k_f·k_e\hfill \\ \hfill \end{array}$$

(3)

where $k_f$ is the fuel consumption per KW·h, and $C_{fuel}$ is the fuel price. In addition to fuel, the operating costs of a ship are expressed as $C_c$. P is the main power of the vessel, and $k_e$ is the CO${}_2$ emission per unit of fuel consumed.

Compared to the raster method [11], based on fields, the environment of an unmanned ship navigating at sea can be portrayed in a more continuous manner. Figure 2 shows the cost potential field; it is on the basis of the map regarding the weather data, according to Equations (1) and (2), and the construction of the cost function is the basis for the navigation module to perform path planning.

3.2. Guidance Field

The guidance field, the design of which is based on the idea of the gravitational field in the APF method [15], is an artificial virtual potential field used to describe the cost of the distance to the end point, and the potential energy value of the gravitational field is specified as a function of being proportional to the distance to the end point to ensure that the potential energy at the end point is the lowest. However, the motion of an unmanned ship has a course, and the steering of the ship will consume a lot of energy, so the concept of direction is added to the design of the guidance field. The guidance field (Figure 3) represents the time penalty of the path point, and this field drawn on the basis of Equation (4):

$$\begin{array}{cc}\hfill & G=\frac{r}{v}·\epsilon ·e^{\frac{{(\psi -\phi )}^2}{{2c}^2}}\hfill \end{array}$$

(4)

where r indicates that the distance between the destination and other points on the map, and v indicates that ship’s speed. $\epsilon$ and c are constants. $\psi$ is the angle between the line connecting the center point of the ship and the target point and the north direction, and $\phi$ is the north angle of the line between each point on the map and the end point. In Figure 3, a small value in the dark area of the color bar represents a low strength of the potential field.

3.3. Field Renewable Local Path Planning Algorithm

3.3.1. Superposition of Fields

From Equations (2) and (3), it can be seen that there were inconsistencies in magnitude between the cost, emission, and the guidance fields, so we needed to normalize the models of all three to solve the problem of the big difference in the scale. So, the normalization treatment was used:

$$\begin{array}{cc}\hfill & \overline{F}=\frac{F-F_L}{F_U-F_L}\hfill \end{array}$$

(5)

In Equation (5), L is lower boundary and U is upper boundary, $\overline{F}$ means the normalized value, $F_L$ means the minimum value number in a series of F values, and $F_U$ means the maximum number in a series of F values. Correspondingly, the value of the potential energy of the superposition field can be calculated by:

$$\begin{array}{cc}\hfill & U=k_1·\overline{EC}+k_2·\overline{EM}+k_3·\overline{G}\hfill \end{array}$$

(6)

In Equation (6), $k_1$, $k_2$, and $k_3$ are the weight coefficients, and $k_3$ is bigger than $k_1$ and $k_2$; EC is the form of fuel cost and charter cost; EM is emission costs; $\overline{EC}$ expresses the EC after normalization; $\overline{EM}$ expresses the EM after normalization; and G expresses the field energy of the guiding field. U expresses the potential energy of the superposition field. The cost field (Figure 2), the emission field (calculated by Equation (3) just like the cost field), and the guidance field are summed by the weights; then, we can obtain the superposition field U, which is shown in Figure 4.

3.3.2. Event Trigger Mechanism

Safety (collision avoidance) is the primary goal of MASS navigation, so we designed an event-triggering mechanism to ensure that MASSs can respond more comfortably when facing collision avoidance situations. The trigger function ET is only related to the distance k from the ship to the obstacle:

$$ET=\left\{\begin{array}{cc}0,\hfill & k>\epsilon \hfill \\ 1,\hfill & k\le \epsilon \hfill \end{array}\right.$$

(7)

As shown in Algorithm 1, steps 12–17, when the radar detects a static or dynamic obstacle within range $\epsilon$, the MASS is transformed from the navigation state to the collision avoidance state, and the event-triggering mechanism gives a signal back to the navigation module to turn off this field, then re-enables it when the ship has passed the obstacle.

3.3.3. Path Planning Based on Superposition Field

In the superimposed field map, the guide field constructs a decreasing gradient to ensure that the field energy at the endpoint is the lowest, while the field energy at the location of the ship is high, and the points in the map at the line of the ship’s location and endpoint are lower than other angle. The cost and emission of the MASS can be expressed as the repulsive potential $P_O(X,Y)$ during the voyage, and the guide field provides the gravitational potential $P_R(X,Y)$.

$$\begin{array}{c}\hfill \left\{\begin{array}{cc}\hfill P_O(X,Y)& =\overline{EC}+\overline{EM}\hfill \\ \hfill P_R(X,Y)& =\overline{G}\hfill \end{array}\right.\end{array}$$

(8)

Superimposing the repulsive potential with the attractive potential, the total potential function is just:

$$P_U(X,Y)=P_O(X,Y)+P_R(X,Y)$$

(9)

The negative gradient of the potential field can be considered as a gravitational force, which guides the ship to move towards the course of the end point. By following the direction of the induced force, the energy of the cost is minimized. Since the ship navigation is attracted by the guiding potential and repelled by the cost potential [30], in a physical sense, total force $\overrightarrow{F_U}$ is a combination of repulsive force $\overrightarrow{F_O}$ and gravitational force $\overrightarrow{F_R}$. After the field strength of the map is determined at a specific moment, the corresponding full gradient direction and contour of the map can be found (Figure 5). Figure 5 shows the gradient map calculated from the superposition field map (Figure 4); in Figure 5, we can observe the gradient direction of each region on the map (indicated by blue arrows). The curves on the map are contour lines, and on the same curve, the potential energy of each point is the same; the darker the color of the curve, the lower its potential energy. The total force on the ship at a determined point (x,y) on the map is:

$$\overrightarrow{F}=-▽P_U(X,Y)$$

(10)

Algorithm 1 is the pseudo-code for local path planning with N steps at a time. In Algorithm 1, we use ▽(p(k)) to denote the gradient direction of the current position and $\alpha$(k) to denote the previous step length of the current cyclic order. In order to make the path obtained from the search conform to the law of ship motion and allow it to be easily transformed into a trajectory, we assume that the value of step size is equal to the value of the ship’s speed $v_s$. The flow of the algorithm can be briefly described as:

• Calculate guide field Gk and initial superposition field $U_k$ according to Equations (1)–(6).
• Judge whether the vessel is in the collision avoidance state. If the distance between the vessel and the obstacle is greater than or equal to danger range $\epsilon$, the vessel is judged to be in the regular navigation state; if the distance between the vessel and the obstacle is less than danger range $\epsilon$, the vessel is judged to be in the regular collision avoidance state. In that case, when the coordinates of the vessel on both the X and Y axes exceed the obstacle or other vessel, the vessel is transferred from the collision avoidance to regular navigation state.
• If the ship is in the state of collision avoidance, calculate gradient ▽(p(k)) from $G_k$; if in the state of regular navigation, calculate the gradient from $U_k$.
• Obtain path points p(k + $\delta$) = p(k) − $\alpha$(k) ∗ ▽(p(k)).

Repeat the above process until the task of path planning is completed and N + 1 path points are obtained.

Algorithm 1 Local path planning.

Input: Gradient of field

Output: Trajectory of MASS p(x)

1: $i=0$   2: $p\left(0\right)=p_{start}$   3: Calculate initial cost and emission field energy after normalization: $E{C}_0$ and $E{M}_0$   4: Calculate initial navigation field energy after normalization: $G_0$   5: Calculate initial superimposition field energy $U_0=k_1·E{C}_0+k_2·E{M}_0+k_3·G_0$   6: Procedure:   7: for$k<N+1$do   8:  Calculate initial cost and emission field energy after normalization: $E{C}_i$ and $E{M}_i$   9:  Calculate initial superimposition field energy $U_i=k_1·E{C}_i+k_2·EMi+k_3·G_i$. 10:  Calculate distance between own ship and obstacle 11:  if distance > $\epsilon$ (Danger range) then 12:   $▽=▽\mathit{U}$; 13:   else, if distance ≤ $\mathbf{\epsilon }$ (danger range) & [$\mathit{x}\left(\mathit{k}\right),\mathit{y}\left(\mathit{k}\right)$] ≤ [$\mathit{o}\mathit{b}\mathit{s}\_\mathit{x}\left(\mathit{k}\right),\mathit{o}\mathit{b}\mathit{s}\_\mathit{y}\left(\mathit{k}\right)$] 14:   $▽=▽\mathit{G}$ 15:   else, if distance ≤ $\mathbf{\epsilon }$ (danger range) & [$\mathit{x}\left(\mathit{k}\right),\mathit{y}\left(\mathit{k}\right)$] > [$\mathit{o}\mathit{b}\mathit{s}\_\mathit{x}\left(\mathit{k}\right),\mathit{o}\mathit{b}\mathit{s}\_\mathit{y}\left(\mathit{k}\right)$] 16:   $▽=▽\mathit{U}$; 17:  end if 18:  while $▽\ne 0$ do 19:   p(k + 1) = p(k) − $\alpha$(k)·▽ (p(k)) 20:  end while 21:  k = k + 1 22:  i = i + 1 23: end for 24: Repeat procedure

4. NMPC-Based Trajectory Tracking Control and Collision Avoidance

In the control module, the accurate trajectory tracking of MASSs needs to be completed; the reference trajectory is transformed from the reference path output by the navigation module in Section 3, which is used as the input of the control module, and the optimal control volumes are output after rolling optimization. In this section, we first impose a time constraint on the input reference path to transform the path point into a trajectory point, and then transform the trajectory into a reference state quantity via Equation (19). We describe the trajectory tracking problem as an error control problem of the ship’s state quantity, design an objective function to minimize the change in state error and control quantity, designing the ship’s collision avoidance behavior as a constraint and embedding the constraint into the MPC controller, and control whether the constraint takes effect by function ET; if ET = 1, this means the constraint is on.

4.1. Maritime Autonomous Surface Ship System Modeling

For the MASS motion model, we refer to the three-degree-of-freedom model of [39]:

$$\begin{array}{cc}\hfill & \dot{x}=ucos\varphi -vsin\varphi \hfill \\ \hfill & \dot{y}=usin\varphi +vcos\varphi \hfill \\ \hfill & \dot{\varphi }=r\hfill \\ \hfill & \dot{u}=\frac{m_{22}}{m_{11}}vr-\frac{d_{11}}{m_{11}}u+\frac{1}{m_{11}}{\tau }_u\hfill \\ \hfill & \dot{v}=-\frac{m_{11}}{m_{22}}ur-\frac{d_{22}}{m_{22}}v\hfill \\ \hfill & \dot{r}=\frac{m_{11}-m_{22}}{m_{33}}uv-\frac{d_{33}}{m_{33}}r+\frac{1}{m_{33}}{\tau }_r\hfill \end{array}$$

(11)

where x and y are the ship’s location, $\varphi$ expresses the ship’s heading angle on the earth-fixed frame, u and v, respectively, express the longitudinal and transverse velocities of the ship’s fixed coordinate system, and r is the angular velocity in yaw around the ship-fixed z-axis. Parameters $m_{11}$, $m_{22}$, and $m_{33}$ express the ship’s Inertia coefficient; $d_{11}$, $d_{22}$, and $d_{33}$ are the hydrodynamic damping coefficients.

We assume that the external environmental forces caused by wind and wave currents are constant with respect to the inertial reference system and that they have three components: surge, sway, and yaw. Therefore, we could modify the ship model written in compact form as [27]:

$$\dot{X}=f\left(X\right)+g\left(X\right)+g_{2b}\left(X\right)w_b$$

(12)

where X = ${\left[xy\varphi uvr\right]}^T$ ∈ ${\mathbb{R}}^6$ is the state vector, u = ${\left[{\tau }_u{\tau }_r\right]}^T$ ∈ ${\mathbb{R}}^2$ expresses the control vector of surge force and yaw moment only, and $w_b$ = ${\left[w_u{w}_v{w}_r\right]}^T$ ∈ ${\mathbb{R}}^3$ expresses the external forces and moment vectors acting in the direction of surge, sway, and yaw and performance in the ship-fixed coordinate system [27]:

$$\begin{array}{c}\hfill f\left(X\right)={\left[\begin{array}{c}ucos\psi -\nu sin\psi \\ usin\psi +vcos\psi \\ r\\ \frac{m_{22}}{m_{11}}vr-\frac{d_{11}}{m_{11}}u\\ n-\frac{m_{11}}{m_{22}}ur-\frac{d_{22}}{m_{22}}\nu \\ \frac{m_{11}-m_{22}}{m_{33}}u\nu -\frac{d_{33}}{m_{33}}r\end{array}\right]}^T\in {\mathbb{R}}^6\end{array}$$

(13)

$$\begin{array}{c}\hfill g\left(X\right)={\left[\begin{array}{cccccc}0& 0& 0& \frac{1}{m_{11}}& 0& 0\\ 0& 0& 0& 0& 0& \frac{1}{m_{33}}\end{array}\right]}^T\in {\mathbb{R}}^6\times {\mathbb{R}}^2\end{array}$$

(14)

$$\begin{array}{c}\hfill g_{2b^{\left(X\right)}}={\left[\begin{array}{cccccc}0& 0& 0& \frac{1}{m_{11}}& 0& 0\\ 0& 0& 0& 0& \frac{1}{m_{22}}& 0\\ 0& 0& 0& 0& 0& \frac{1}{m_{33}}\end{array}\right]}^T\in {\mathbb{R}}^6\times {\mathbb{R}}^3\end{array}$$

(15)

By using a virtual vessel that has the same dynamics described as Equation (11), we generated time-varying predicted trajectories, but the virtual vessel was independent of external forces:

$$\begin{array}{cc}\hfill & \dot{X_r}=f\left(X_r\right)+g\left(X_r\right)u\hfill \end{array}$$

(16)

The time series of state and control vectors of MASS are defined as follows:

$$\begin{array}{cc}\hfill & X\left(t\right)\in X\hfill \\ \hfill & u\left(t\right)\in U\hfill \end{array}$$

(17)

Here, X and U are sets of state and control vectors.

Combining Equations (16) and (17), we can obtain the time-varying prediction equation for the MASS state variables:

$$\begin{array}{c}\hfill X\left(t+1\right)=f\left(X\left(t\right),u\left(t\right)\right),t>0\end{array}$$

(18)

where $x\left(t\right)$ ∈ ${\mathbb{R}}^6$ expresses the vector of state variables; $u\left(t\right)$$\in {\mathbb{R}}^2$ expresses the vector of control variables.

4.2. Trajectory Tracking Control

The original sequence of the input signal in this section is the reference path point $p\left(x\right)$ obtained in Section 3.3, and by the time constraint, we can transform it into reference trajectory point $p\left(t\right)$. Consider the reference trajectory $p\left(t\right)=[x_d\left(t\right),y_d\left(t\right)]$, which defines the position of MASSs on the local horizontal plane. Here, we assume that the following conditions hold [32]:

Assumption 1: Trajectory $p\left(t\right)$ and its derivatives are expected to be smooth and bounded, satisfying: $0\le p\le {∥p\left(t\right)∥}_{\infty }\le \stackrel{-}{p},0\le \underline{p_1}\le {∥\dot{p\left(t\right)}∥}_{\infty }\le \stackrel{-}{p_1}<\infty ,0\le \underline{p_2}\le {∥\ddot{p\left(t\right)}∥}_{\infty }\le \stackrel{-}{p_2}<\infty ,0\le \underline{p_3}\le {∥\stackrel{⃛}{p\left(t\right)}∥}_{\infty }\le \stackrel{-}{p_3}<\infty .$

To correspond to the six state quantities of MASSs in Equation (11), we can extend the desired trajectory $p\left(t\right)$ to a reference state quantity $X_d\left(t\right)$ by Equation (17), where $X_d\left(t\right)=[x_d\left(t\right),y_d\left(t\right),\varphi {}_d\left(t\right),u_d\left(t\right),v_d\left(t\right),r_d\left(t\right)]$, and the four extended reference states are calculated by the following equations:

$$\begin{array}{cc}\hfill & {\varphi }_d\left(t\right)=atan2\left({\dot{y}}_d\left(t\right),{\dot{x}}_d\left(t\right)\right)\hfill \\ \hfill & u_d\left(t\right)=\sqrt{{\dot{x}}_d^2\left(t\right)+{\dot{y}}_d^2\left(t\right)}\hfill \\ \hfill & v_d\left(t\right)=0\hfill \\ \hfill & \left.r_d\left(t\right)=\left({\dot{x}}_d\left(t\right){\ddot{y}}_d\left(t\right)-{\dot{y}}_d\left(t\right){\ddot{x}}_d\left(t\right)\right)/{\dot{\left(x\right.}}_d^2\left(t\right)+{\dot{y}}_d^2\left(t\right)\right)\hfill \end{array}$$

(19)

where $atan2$ is the four-quadrant arctangent operator in Matlab, and $\dot{f(·)}$ denotes the derivative. Then, it can be verified that $X_d\left(t\right)$ satisfies kinematic Equation (11).

In this way, we obtain the expressions of all reference states of MASSs. In the control module, the nonlinear MPC framework uses a time-varying motion model running online to predict the vessel’s motion state. Solving a finite-time domain nonlinear optimal control problem can be regarded as an optimization problem, and the prerequisite for solving optimization problems is to establish a suitable optimization objective function. In order to make MASSs fast-track the expected trajectory, we use the state error and control error of the finite-time domain to determine the stage cost; in addition, the vessel can contribute a terminal cost function $X\left(T\right)$ to ensure the stability of the control system [29]. Then, the expression of the objective function is given by:

$$\begin{array}{cc}\hfill & \underset{\widehat{u}\in S\left(\delta \right)}{min}J\left(X,u\right)={\int }_0^T{∥\stackrel{\sim }{X}\left(s\right)∥}_Q^2+{∥\stackrel{\sim }{u}\left(s\right)∥}_R^{2}s+{∥\stackrel{\sim }{X}\left(T\right)∥}_p^2\hfill \end{array}$$

(20)

where $\widehat{X}\left(s\right)$ is the predicted ship trajectory based on control variables $\widehat{u}\left(s\right)$, which evolves from $x\left(t_0\right)$ using the ship motion model; $\stackrel{\sim }{\mathbf{X}}=\widehat{\mathbf{X}}-{\mathbf{X}}_{\mathbf{d}}$ expresses the state variable error and $\stackrel{\sim }{\mathbf{u}}=\widehat{\mathbf{u}}-{\mathbf{u}}_{\mathbf{d}}$ is the control variable error; $S\left(\delta \right)$ expresses the group of segmented constant functions characterized by the sampling period $\delta$; and T = N $\delta$ indicates the control interval. Q, R, and P express the weight matrices, and all are positive definite matrices.

Based on the dynamic principle of MASSs, we have the following constraint:

$$\begin{array}{cc}\hfill s·t·& \dot{\widehat{x}}=f\left(\widehat{x}\left(s\right),\widehat{u}\left(s\right)\right)\hfill \\ \hfill & \widehat{x}\left(0\right)=x\left(t_0\right)\hfill \\ \hfill & u_{min}<\widehat{u}\left(s\right)<u_{max}\hfill \end{array}$$

(21)

We treat collision avoidance as a controller-dependent planning problem, but due to unfavorable environmental conditions and neglected ship dynamics, the controller may not be able to achieve a solution to the problem. To overcome this, we utilize the ability of NMPC to systematically integrate state constraints into its optimization problem, and assuming that the obstacles are circular, autonomous collision avoidance can be transformed into an inequality constraint integrated into the design of the controller with the following [29]:

$$\begin{array}{cc}\hfill & sqrt\left({\left(x\left(t\right)-x_0^i\left(t\right)\right)}^2+{\left(y\left(t\right)-y_0^i\left(t\right)\right)}^2\right)\ge R+R_0^i\hfill \end{array}$$

(22)

where $x\left(t\right)$ and $y\left(t\right)$ are the predicted positions of one’s own vessel over the prediction horizon, $x_0^i\left(t\right)$ and $y_0^i\left(t\right)$ are the position of ith dynamic or statistic obstacles, and R and $R^i$ are the hazard radii of one’s own vessel and obstacle.

Eventually, the nonlinear optimization control problem of MASS trajectory tracking can be integrated as:

$$\begin{array}{cc}\hfill & \underset{\widehat{u}\in S\left(\delta \right)}{min}J\left(x,u\right)={\int }_0^T{∥\stackrel{\sim }{x}\left(s\right)∥}_Q^2+{∥\stackrel{\sim }{u}\left(s\right)∥}_R^{2}s+{∥\stackrel{\sim }{x}\left(T\right)∥}_p^2\hfill \\ \hfill s·t·& \dot{\widehat{x}}=f\left(\widehat{x}\left(s\right),\widehat{u}\left(s\right)\right)\hfill \\ \hfill & \widehat{x}\left(0\right)=x\left(t_0\right)\hfill \\ \hfill & u_{min}<\widehat{u}\left(s\right)<u_{max}\hfill \\ \hfill & sqrt\left({\left(x\left(t\right)-x_0^i\left(t\right)\right)}^2+{\left(y\left(t\right)-y_0^i\left(t\right)\right)}^2\right)\ge R+R_0^i\hfill \end{array}$$

(23)

The flow of the NMPC-based MASS trajectory tracking algorithm can be briefly described as follows:

• The optimization problem (23) is solved with the current state error $\tilde{X}\left(t_0\right)$ as the initial state and $u(t_0+k\delta )=\left[{\tau }_u(t_0+k\delta ),{\tau }_r(t_0+k\delta )\right]$, (k = 0, 1, 2,…, N − 1) as a set of optimal solutions to the control problem.
• The system only uses the first optimal solution as the control output for one sampling period $[t_0,t_0+\delta ]:u\left(t_0\right)=[{\tau }_u\left(t_0\right),{\tau }_r\left(t_0\right)]$, as a control decision for the MASS.
• At moment $t_0+\varphi$ at the beginning of the next sampling period, the optimization problem (23) is solved again using the latest measurement data $\tilde{X}(t_0+\delta )$.

The above process is repeated until the MASS tracking control task is completed.

5. Collaborative Optimization of Path Planning and Tracking Control

In the traditional research of path planning and trajectory tracking control, the information transfer between navigation and control modules is often offline and one-way, which means it is non-closed-loop information transfer. This transfer mode cannot achieve online planning and control, and for the ship motion system with large inertia and time lag, online planning is more effective than offline planning. Therefore, we designed a collaborative optimization framework of path planning and trajectory tracking control to complete online planning and control and achieve the purpose of the collaborative optimization of navigation and control modules.

Based on the characteristics of NMPC data sampling, the tracking control can be easily combined with the path planning method of rolling planning proposed in the previous section, and the information transfer between the path planning module and the control module can be easily accomplished at each sampling moment through the information transfer interface. The closed-loop structure of the autonomous navigation system, composed of path planning and NMPC trajectory tracking control modules, is depicted in Figure 6.

5.1. Interface I for Path Planning

The output of the NMPC optimization framework is the control variable $u=[{\tau }_u,{\tau }_r]$, and according to the time-varying prediction model of the MASS state volume, we can calculate $X_R(t+\delta )=f(X_R\left(t\right),u_R\left(t\right))$ and matrix $I\in {\mathbb{R}}^6$ as the interface for information transfer between the control and navigation modules to obtain the initial parameters of the path planning algorithm:

$$\begin{array}{cc}\hfill & S_I\left(0\right)=I\times X_R(t_0+\delta )\hfill \end{array}$$

(24)

In Equation (24), $I(1,1)=1,I(2,2)=1,I(3,3)=1$, and the other values in the matrix = 0.

5.2. Process of the Closed-Loop Structure

The navigation module can generate N + 1 geometric path points after completing a path planning task at time t. Due to the design feature of the search step being equal to the ship’s speed (Section 3.3), it is easy to convert the generated geometric path into a time-dependent trajectory profile. The reference trajectory can be expressed as:

$$\begin{array}{cc}\hfill & \left\{p\right(x),x\in [t,t+(N+1)·\delta ],xisanintegermultipleof\delta \}\hfill \end{array}$$

(25)

$p\left(x\right)$ can be converted to the reference state quantity in real time via Equation (19), just as follows:

$$\begin{array}{cc}\hfill & X_d={\left[x_d,y_d,{\varphi }_d,u_d,{\nu }_d,r_d\right]}^T\hfill \end{array}$$

(26)

Then, state error $\tilde{X}\left(t_0\right)$ is used as the input of the control module, and after the optimal solution of NMPC, the optimal control output sequence $U\left(T\right),T\in (t_0,t_0+N\delta )$ is obtained, and the first control variable $u\left(t_0\right)={[{\tau }_u\left(t_0\right),{\tau }_r\left(t_0\right)]}^T$ is taken as the control command of MASS at moment $t_0$. According to Equation (24), the partially predicted state at moment $t_0+\delta$ can be calculated by optimal control $\tau =u\left(t\right)$:

$$\begin{array}{cc}\hfill & X_R\left(t_0+\delta \right)=f\left(x_R\left(t_0\right),u_R\left(t_0\right)\right)\hfill \end{array}$$

(27)

$$\begin{array}{cc}\hfill & S_I=\left[x\left(t_0+\delta \right),y\left(t_0+\delta \right),u\left(t_0+\delta \right)\right]=I\ast X_R\left(t_0+\delta \right)\hfill \end{array}$$

(28)

$S_I$ is used as input to the navigation module for the path search for the next sampling period $[t+\delta ,t+2\delta ]$, and the potential field is regenerated based on the updated starting point and the ship’s speed:

$$\begin{array}{cc}\hfill & x\left(t+\delta \right)=\left[x\left(t+\delta \right),y\left(t+\delta \right)\right]\hfill \\ \hfill & U=U\left(u\left(y+\delta \right),\dots \right)\hfill \\ \hfill & G=G\left(x\left(t+\delta \right),y\left(t+\delta \right)\dots \right)\hfill \end{array}$$

(29)

In the updated potential field environment, path planning is performed for $t+\delta$ moments to form a new series of reference trajectories. The above process is repeated until the autonomous navigation task of the MASS is completed. The rolling-optimization-based path planning and NMPC-based trajectory tracking combined with dynamic and static collision avoidance control algorithms are summarized in Algorithm 2:

Algorithm 2 Path planning and trajectory tracking control.

Input: Field and reference trajectory

Output: Optimum control and real trajectory

1: $i=0$   2: initial field $U_0$   3: initial gradient matrix $G_0$   4: Procedure   5: for$i=1:sim\_tim$do   6:  for $m=1:N+1$ do   7:   use path planning Algorithm 1 to calculate N + 1 reference position $p\{N+1\}$   8:  end for   9:  for $k=1:N$ do 10:   Transform reference position to reference trajectory and reference state $X_d\left(i\right)$ 11:   Use NMPC module to calculate control quantity $\tau$ 12:  end for 13:  Calculate next state $X_R(i+1)$ by $\tau$; $S_I(i+1)=I\times X_R(i+1)$; 14:  initial field $U=U_i$, initial gradient matrix $G=G_i$; 15:  $i=i+1$; 16: end for 17: Repeat procedure

6. Simulation

6.1. Experimental Design

6.1.1. Marine Environment

Marine meteorological data can be obtained through many channels, and in this paper, in order to simplify the process of data acquisition, multi-point fitting was used to acquire more meteorological data. In the actual marine meteorological environment, wind and wave data can be obtained separately. Due to the different timeliness of data, wind and wave data cannot be obtained accurately at the same time; to ensure the accuracy of wave data, wave height is calculated from wind speed information using the following equation [12]:

$$\begin{array}{cc}\hfill & h=\frac{0.7·{\left(\frac{gF}{{\nu }_{wind}^{2}d}\right)}^{\frac{1}{3}}·{\nu }_{wind}^2}{g}\hfill \end{array}$$

(30)

In Equation (30), gravity acceleration $g=9.8$ m/s${}^2$, and ${\nu }_{wind}$ is wind speed. F is the length of wind area, which refers to the region of sea where the wind strength conditions are almost the same.

Since the magnitude of ocean currents is less affected by other factors, the current speed is simply set to a constant value throughout the map. The final constructed environment for navigation is shown in Figure 7. Figure 7 shows the numerical overlay of wind speed, wave height, and current velocity, where wind speed data are obtained via cubic spline interpolation fitting, wave height data are calculated via Equation (30), and current velocity is a constant value. The brighter the color of the area (the higher the value), the more severe the weather environment.

6.1.2. Experimental Parameters

In the simulation experiments, we chose a map range of 100 nm × 100 nm and a KCS container ship. The hydrodynamic parameters and other ship parameters are shown in

Table 2. KCS ship and weight parameters.

Parameters	Value	Unit
$m_{11}$	897.77	kg
$m_{22}$	1718	kg
$m_{33}$	105.34	kg
$d_{11}$	42.77	kg/s
$d_{22}$	863	kg/s
$d_{33}$	51.88	kg/s
S (wetted  surface)	9424	m${}^2$
A (wind  area)	2200.43	m${}^2$
L	230	m
B	32.2	m
$k_1$	1
$k_2$	1
$k_3$	16
$C_w$	0.57e${}^{-3}$
$\nu$	1.27e${}^{-6}$	m${}^2$/s
$C_a$	0.1

. A predicted horizon length NP of 420 min and a sampling interval of 4 min were chosen in order to reduce the computation, so the prediction step Np was 105 and the control step Nc was 15.

The controller weights were as follows: Q(1,1) = 1 × 10${}^1$; Q(2,2) = 1 × 10${}^1$; Q(3,3) = 10; Q(4,4) = 1; Q(5,5) = 1; Q(6,6) = 1; R(1,1) = 1 × 10${}^{-2}$; R(2,2) = 1 × 10${}^{-2}$; R1(1,1) = 1 × 10${}^{-2}$; R1(2,2) = 1 × 10${}^{-2}$.

The vessel control constraints of the autonomous navigation and control system are defined as follows (N·m):

$$\left[\begin{array}{c}-2000\\ -2000\end{array}\right]\le U\le \left[\begin{array}{c}2000\\ 2000\end{array}\right]$$

(31)

For simplicity of calculation, the value of drag coefficient $C_{ts}$ was made to approximate frictional resistance $C_F$, without considering the residual resistance. The calculation equation is as follows [40]:

$$C_F=\frac{0.075}{{\left(log\left(Re\right)-2\right)}^2},Re=\frac{v_s·L}{\nu }$$

(32)

where $v_s$ is the speed of vessel over water, L expresses length of the vessel, and $\nu$ expresses the viscosity of water.

6.1.3. Experimental Design

In this paper, simulation experiments were performed in MATLAB, and four sets of experiments were conducted.

Experiment 1: A comparison experiment between single-objective planning with cost reduction and multi-objective planning with cost and emission reduction was designed to show the impact of adding carbon emissions to path planning.

Experiment 2: A comparison experiment was conducted between the RHNC and Offline-SPNC systems under conventional conditions to validate the effectiveness of the RHNC system in reducing economic costs and emissions.

Experiment 3: Autonomous collision avoidance experiments were conducted with the MASS controller in static collision avoidance and pair encounter scenarios to verify the control and collision avoidance performance of NMPC.

Experiment 4: A MASS autonomous navigation experiment containing collision avoidance state and regular state transition was conducted to verify the effectiveness of the event-triggering mechanism.

6.2. Simulation Results

6.2.1. Contrast Experiments Considering Carbon Emissions and No Carbon Emissions in Local Path Planning

We first conducted a series of controlled experiments to verify the effect of the local path planning algorithm presented in the previous article after adding carbon emissions as a penalty function. This set of experiments was divided into two parts: path planning considering economic cost and carbon emission, and path planning considering only economic cost without carbon emission. The coordinates of the starting point of both parts of the experiment were O(0, 0), the coordinates of the local target point were D(100, 100), the search step length of the path planning was set to 1.2 nautical miles, and the total step length NL was 110. Figure 8 shows that the paths searched by the path planning algorithm offline (blue paths) become more tortuous than the paths searched by the algorithm considering only the economic cost objective (red paths) after introducing the carbon emission cost. We can see that the potential energy of the path points obtained by the search algorithm was lower after the emission field in Figure 9 was added, which reflects the effect of emission field construction and field superposition. This shows the effectiveness of our proposed route planning algorithm in considering carbon emissions, and that the navigation module takes into account the goal of reducing costs and emissions when planning paths.

6.2.2. Simulation Experiments of Receding Horizon Navigation and Control System in Global Planning

In this subsection, we describe the simulation experiments conducted to verify the effectiveness of the proposed receding horizon navigation and control (RHNC) system for global navigation and the control of MASSs in Section 5. In the experiments, the step size of the path planning algorithm in the navigation module was 16, and the control step size of the MPC in the control module was 15. Figure 10 shows the map field strengths at six sampling moments with sim-time (Np) of 20, 40, 60, 80, 100, and 120 and the paths that MASSs navigated at these six sampling moments. We can find from the graph that the map field intensity changed at different sampling moments, which was caused by changes in ship speed and the angle between the ship’s location and the target point, which verifies the changes caused by the rolling optimization framework to the path planning algorithm. After the simulation experiment, we can verify that the proposed navigation and control modules can plan the navigation path according to the field strength on the map while controlling the autonomous ship to follow the local path planned by the navigation module, and the navigation module can obtain the state information of the autonomous ship in real time to plan the local path in the next sampling moment.

Figure 11 shows a time-varying diagram of control quantity $\tau$. From the diagram, we can see that ${\tau }_u$ and ${\tau }_r$ fluctuated greatly in the range of 240–300 min, and then returned to normal, which means that the navigation module of the autonomous vessel plans a route with a relatively large turning angle during this time, and the control module controls the autonomous vessel to make a larger turn after obtaining the reference trajectory input. This phenomenon is verified in Figure 10.

In order to further verify the effectiveness of the navigation and control system with rolling optimization proposed in this paper, we designed a control experiment using offline shortest path planning with NMPC trajectory tracking control (Offline-SPNC); a trajectory diagram is shown in Figure 12.

For simplicity of calculation, we took the average of the magnitude of the economic cost field energy and emission field energy (before normalization treatment) of the path point at two consecutive sampling moments as the increase rate of economic costs and emissions for one consecutive time period. The economic costs and emissions of the navigation and control systems are shown in

Table 3. Costs and emissions.

Navigation and Control System	EC ${10}^4$ (USD)	EM (t)	EM per n mile (USD)	EM per n mile (kg)	Time (min)
RHNC	1.047	31.2642	80.74	241.10	480
Offline-SPNC	1.1656	37.8385	88.81	288.31	460

. The approximate coordinates of the final position A1 of a MASS under the control of the Offline-SPNC system are approximately (92.8 n mile, 92.8 n mile), and those of the final position A2 of a MASS at the end of the simulation of the RHNC system are (94.3 n mile, 89.8 n mile), and the linear mileage of both is approximately 131.24 n mile for Offline-SPNC and 129.67 n mile for RHNC. Compared with the Offline-SPNC system, per n mile of the RHNC system, the economic cost is reduced by about 9.08%, and the emissions are reduced by about 16.37%.

6.2.3. Simulation Experiments of Real-Time Path Planning and Tracking Control with Collision Avoidance

To evaluate the performance of the NMPC controller used to avoid collisions between static obstacles and dynamic objects, the control performance of one static scenario and one dynamic scenario was analyzed. The predicted trajectory of the target ship was calculated using an isokinetic model under the assumption that the other ship does not have a collision avoidance scenario. Finally, a comprehensive simulation experiment containing all elements (navigation and control, static collision avoidance, dynamic collision avoidance, and event-triggering mechanism) was performed.

Scenario: Static Collision Avoidance

The map range used in this group of experiments was 10 n miles × 10 n miles, the radius of the static obstacle hazardous collision avoidance range was set to 2 n miles, the chosen sampling interval was 20 s in consideration of reducing computation, the total prediction step $Np$ was 110, that is, the total simulation time was 2200 s, and the control step of NMPC controller $N_c$ was 15.

The change process of the ship’s trajectory and control variables are shown in Figure 13 and Figure 14.

Scenario: Head-on

In the experiments of dynamic collision avoidance, the dangerous collision avoidance radius of the ship was set to 1 n mile, and the chosen sampling interval was 20 s in consideration of reducing computation. The total prediction step Np in this scenario was 80, so that the total simulation time was 1600 s, and the control step of the NMPC controller Nc was 15.

Supposing a ship is located at $(X_T,Y_T)$, its trajectory is as follows:

$$\left\{\begin{array}{c}{x}_T\left(t\right)=4-0.05·t\\ y_T\left(t\right)=4-0.05·t\end{array}\right.$$

(33)

where t represents the frequency of the sampling moment rather than the actual time T, that is, $t=T/\delta$. For each sampling moment interval of 20 s, Equation (32) means that every 20 s along the X and Y axes, the obstacle ship travels a distance of 0.5 n miles.

The trajectory diagram is shown in Figure 15, and the surge force and yaw moments are expressed in Figure 16. In this case, the collision avoidance constraint and the longitudinal swing force and transverse swing moment constraint are satisfied, and after the collision avoidance is completed, the rolling optimized navigation module will plan a new path according to the state of the intelligent merchant vessel; then, the NMPC controller in the module will control the vessel to track the new trajectory. Compared with offline planning, a MASS under the RHNC system can save part of the journey to return to the original path.

6.2.4. RHNC System with Collision Avoidance and Event Triggering

In order to verify the effectiveness of the whole system, a comprehensive simulation experiment was conducted. Supposing a ship is located at $(X_T,Y_T)$, its trajectory is as follows:

$$\left\{\begin{array}{c}{x}_T\left(t\right)=50-0.75·t\\ y_T\left(t\right)=40-0.75·t\end{array}\right.$$

(34)

where t represents the frequency of the sampling moment rather than the actual time T, that is, $t=T/\delta$. For each sampling moment interval of 4 minutes, Equation (33) means that every 4 minutes along the X and Y axes, the obstacle ship travels a distance of 0.75 n miles.

Assuming an island is located at (75 n mile, 65 n mile) with an impact radius $R_{o1}$ of 3 n miles, the hazard radius $R_s$ of the ship is set to 1.5 n miles, and the hazard warning distance of the event trigger mechanism ${\epsilon }_1$ = 4.5 n mile (distance between a MASS and a dynamic obstacle) and ${\epsilon }_2$ = 6 n mile (distance between a MASS and a statistic obstacle).

In the experiment, the step size of the path planning algorithm in navigation module $N_L$ was 16, and the control step size of MPC in control module $N_C$ was 15. Figure 17 shows the field strength of the six sampling moments with a simulation step ($Np$) of 20, 40, 60, 80, 100, and 120 and the trajectories that a MASS sailed at these six sampling moments. In simulation steps 1–20, the MASS sailed in the conventional state considering “cost and emission”; in simulation steps 20–40, the ship found the other ship and changed from the conventional state to the collision avoidance state, and then back to the conventional state again after the collision avoidance was completed. In simulation steps 40–60 and 60–80, autonomous navigation continued in the conventional state considering “cost and emissions”. In simulation steps 80–100 and 100–120, the ship found the static obstacle and changed from the regular state to the collision avoidance state, and then back to the regular state after avoiding and crossing the obstacle (the trajectory of the obstacle beyond the map range was not visualized). After verifying the simulation experiment, the navigation and control modules proposed in this paper can cooperate to optimize the system for real-time planning and control. In the process of control, real-time collision avoidance can be accomplished, and the event trigger mechanism can issue safety warnings for the RHNC system to keep the ship safe in the first place. Figure 18 shows a time-varying diagram of control quantity $\tau$ in the process.

7. Conclusions

In this paper, a new navigation control scheme was designed to meet economic and green objectives in navigation and safe trajectory tracking in control. The path search of the navigation module was based on the path planning method of potential field theory, and the idea of potential field construction was derived from APF [15]. The guidance field and cost emission field were constructed based on the gravitational field and repulsive field, respectively, in APF and were integrated with NMPC with rolling optimization characteristics. In addition, the closed-loop structure of information transfer was studied to design a rolling optimization navigation and control system framework. From the experimental results, it can be seen that the RHNC system proposed in this paper can plan energy-saving and emission reduction paths for MASSs in wind and wave environments and can control an intelligent merchant ship to follow the planned path. Meanwhile, through the design of an event-triggering mechanism, the system can close the cost emission field and change into a collision avoidance state when detecting obstacles or other ships. The experimental results of static and dynamic collision avoidance scenarios can also prove that the system has good dynamic and static collision avoidance performance and can re-route for MASSs after confirming that the ship is in a safe state. In our research, all obstacles are circular, circular domains are quite old. Therefore, there are some potential limitations of this paper. In future research, we will use newer and more scientific ship domains models to design the MPC controller’s collision avoidance constraints.