# Parameters’ Identification of Vessel Based on Ant Colony Optimization Algorithm

Authors: Chen Zhao; Xiaojian Li
Journal: Mathematical Problems in Engineering
 (2021)

---

## Abstract

In this paper, the ant colony optimization (ACO) method is used to identify the parameters of a 3-DOF nonlinear vessel model. Identifying the parameters is abstracted as a nonlinear optimization problem to solve through the ant colony optimization algorithm. The identification procedure is divided into two parts. The first part of the identification procedure is to identify the parameters related to surge motion. The second part of the identification procedure is to identify the rest parameters of the vessel’s kinetics model. In the surge model identification procedure, the transient motor speed is used to generate the training data, and in the sway and yaw motion identification procedure, the zigzag maneuvering with different motor speeds is used to generate the training data. All the parameters are identified by the ACO method and the least-square (LS) method based on the training data and then validated on the validation data. The prediction performance of parameters identified by different methods is compared in the simulation to demonstrate the effectiveness of the ACO algorithm.


---

## Body

## 1. Introduction

Unmanned technology has developed rapidly in recent years and attracted more and more attention from academia and industry. Many kinds of unmanned products such as unmanned vehicles, unmanned aerial vehicles, unmanned surface vessel, and unmanned submersibles have been widely deployed as a network in various situations such as scientific research, environmental missions, ocean resource exploration, military use, and other applications [1] to help people improve work efficiency. As an intelligent device that works on the water, the unmanned surface vessel can work as a node of the overall unmanned network and extend the working range of the entire unmanned network to the surface of the water and underwater. The level of intelligence of the unmanned network is related to the level of intelligence of every node in the network. A sufficiently intelligent unmanned surface vessel contains a lot of techniques such as navigation techniques [2, 3], guidance techniques [4, 5], and control techniques [6–8].The modeling of the unmanned surface vessel is imperative for both control method design and simulation study purposes [1]. To describe the vessel motion in the surge, sway, and heave, a nonlinear 6-DOF model with different dimension motions coupled together is imperative. However, the nonlinear 6-DOF model has lots of parameters which are very difficult to identify, the vessel algorithms designed based on the nonlinear model are very difficult, and the real-time performance of the algorithms is hard to guarantee. For algorithm design, some simple models such as the first-order Nomoto model proposed by Nomoto in 1957 [9, 10] are used. The Nomoto model describes vessel motion approximately under the assumption that the forward speed of the vessel changes slowly. Except for the Nomoto model, the 3-DOF model, including surge velocity, sway velocity, and yaw angular velocity, is also widely used in vessel control [11–13].The prediction results of the Nomoto model and 3-DOF model cannot meet the requirements of the simulation study. How to identify the parameters of the nonlinear 6-DOF model is a challenging problem. The towing carriage can identify the parameters of the vessel model such as CyberShip II [14], but the size of the towing carriage limits the size of the identified vessel. Some research use CFD [15] to simulate the ship motion and to identify the parameters [16], but the accuracy of the CFD in a complex situation is still worthy of further research. Many research studies choose full-scale vessel trials to identify the vessel model; the approach to model identification involved adapting model parameter value which was proposed in [17]. In [18], a new transformed multi-innovation least-squares (TMILS) algorithm is developed; the model structure is identified first and then the parameters are identified with full-scale trial data. An artificial neural network was used in [19] to model a high-speed craft with sea trial data.As for the identify method, the least-squares method [20–23] has been a widely used method in system identification procedure, but the least-squares method can only deal with the linear problem; some identified procedure are abstracted as a nonlinear optimization problem. Support vector machines have the advantage in solving the nonlinear optimization problem, and it was used to identify the parameters of the vessel in [24–27]. Except for the least-square method and support vector machines method, some other methods are also used to identify the parameters. In [28], the identification procedure was divided into two steps: the first step was to determine the structure of the nonlinear model and the second step was parameter estimation refinement by using a nonlinear prediction error method with the unscented Kalmen filter. In [29], a sensitivity analysis and SQP method were used to identify the hydrodynamic coefficients of the Esso Bernicia tanker.In this paper, the ant colony optimization method is used to identify the parameters of the vessel kinetics model. The parameter-identified problem is summarized as a nonlinear least-square problem, and the ant colony optimization method is used to solve the nonlinear least-square problem. The solution to the problem is the parameters to be identified. The identification procedure is divided into two parts [23]: the parameters related to surge motion were identified first by an experiment and then the rest of the parameters related to sway and yaw motion were identified by the zigzag test.The organization of this paper is as follows. In Section2, a nonlinear discrete 3-DOF state-space model of USV with surge speed, lateral speed, and yaw angle as state variables is established. Section 3 discusses the ant colony optimization method. Section 4 illustrates and analyzes the identification results. Section 5 summarizes the conclusions.


## 2. Problem Formulation

### 2.1. 3-DOF Nonlinear Model

In this section, the 3-DOF nonlinear symmetric model is established. By neglecting the heave, roll, and pitch motion, we consider the maneuvering models proposed by [17, 30](1)Mν˙+Cνν+Dνν=τ+τwind+τwave.(i)
Vectorν=uvr⊤ denotes vessel surge velocity, sway velocity, and yaw angular velocity, respectively.(ii)
MatrixM accounts for inertial effects:(2)M=m−Xu˙000m−Yv˙mXg−Yr˙0mXg−Yr˙Iz−Nr˙.(iii)
MatrixCν accounts for centrifugal and Coriolis effects:(3)Cν=00c1300c23c13−c230,where(4)c13=−mXgr+v+Yv˙v+12Nv˙+Yr˙r,c23=mu−Xu˙u.(iv)
MatrixDν accounts for viscous and dissipative effects:(5)Dν=d11000d22d230d32d33,where(6)d11=−Xu−Xuuu−Xuuuu2,d22=−Yv−Yvvv−Yrvr,d23=−Yr−Yvrv−Yrrr,d32=−Nv−Nvvv−Nrrr,d33=−Nr−Nvrv−Nrrr.(v)
Vectorτ denotes the force and torques generated by actuators.(vi)
τwind denotes force and torques caused by wind(vii)
τwave denotes force and torques caused by wave.The vessel model used in this paper is CyberShip II [14], which has a bow thruster, two thrusters, and two rudders. For the convenience of modeling force and torques, the bow thruster is excluded from the actuator model, and two thrusters and rudders are equivalent to one. Thus, the actuator model can be written as(7)τ=Bτactν,n,δ,where matrix B is the actuator configuration matrix and vector τactν,n,δ is the force and torques related to speed of the motor n and rudder angle δ:(8)Bν=20020−2lxR,τactν,n,δ=TL⊤,with(9)T=Tnnnn−Tnunu,L=Lδδ−Lδδδδuu.


### 2.2. Discrete Nonlinear Model

In the vessel experiments, the experiment data are collected in the discrete form. Also, vessel motion controllers are usually designed in the discrete zone. Therefore, the continuous-time model equation (1) is approximated as a discrete model, and the parameters of the discrete model are identified in this paper. The external force is all ignored for the convenience of discretizing the continuous-time model, and equation (1) can be rewritten as(10)Mν˙=Acα,where the matrix Ac and vector α denote the time-varying terms and constant parameters in the vessel model:(11)Acα=−Cνν−Dνν+Bτactν,n,δ,α=uuuu3r2vrn2nu000uvurvrvvvrrvrruuδuuδδuvurvrvvvrrvrruuδuuδδ⊤,Ac=XuXuuXuuumxg−12Nv˙+Yr˙m−Yv˙2Tnn−2Tnu0000Xu˙−mYvYr|Yvv|Yvr|YrvYrr2Lδ−2LδδYv˙−Xu˙12Nv˙+Yr˙−mxgNvNrNvvNvrNrvNrr−2lxRLδ2lxRLδδ⊤.The sample time of collecting data in simulation is denoted asTs; equation (10) can be approximated by the back Euler integration method at time k as(12)νk+1=νk+TsM−1Acαk=Bcαk+TsM−1Acαk,where Bc=100000000000100000000001000000⊤.Because the matrixBc and TsM−1Ac are both constant matrix, thus, equation (12) can be rewritten as(13)uk+1=a1uk+a2ukuk+a3uk3+a4rk2+a5vkrk+a6nk2+a7nkuk,(14)vk+1=b1ukvk+b2ukrk+b3vk+b4rk+b5vkvk+b6vkrk+b7rkvk+b8rkvk+b9ukukδk+b10ukukδkδk,(15)rk+1=c1ukvk+c2ukrk+c3vk+c4rk+c5vkvk+c6vkrk+c7rkvk+c8rkvk+c9ukukδk+c10ukukδkδk.The parametersa1,…,a7, b1,…,b10, and c1,…,c10 are unknown parameters. In Section 3, an ACO method is used to identify parameters of equations (13)–(15).


## 2.1. 3-DOF Nonlinear Model

In this section, the 3-DOF nonlinear symmetric model is established. By neglecting the heave, roll, and pitch motion, we consider the maneuvering models proposed by [17, 30](1)Mν˙+Cνν+Dνν=τ+τwind+τwave.(i)
Vectorν=uvr⊤ denotes vessel surge velocity, sway velocity, and yaw angular velocity, respectively.(ii)
MatrixM accounts for inertial effects:(2)M=m−Xu˙000m−Yv˙mXg−Yr˙0mXg−Yr˙Iz−Nr˙.(iii)
MatrixCν accounts for centrifugal and Coriolis effects:(3)Cν=00c1300c23c13−c230,where(4)c13=−mXgr+v+Yv˙v+12Nv˙+Yr˙r,c23=mu−Xu˙u.(iv)
MatrixDν accounts for viscous and dissipative effects:(5)Dν=d11000d22d230d32d33,where(6)d11=−Xu−Xuuu−Xuuuu2,d22=−Yv−Yvvv−Yrvr,d23=−Yr−Yvrv−Yrrr,d32=−Nv−Nvvv−Nrrr,d33=−Nr−Nvrv−Nrrr.(v)
Vectorτ denotes the force and torques generated by actuators.(vi)
τwind denotes force and torques caused by wind(vii)
τwave denotes force and torques caused by wave.The vessel model used in this paper is CyberShip II [14], which has a bow thruster, two thrusters, and two rudders. For the convenience of modeling force and torques, the bow thruster is excluded from the actuator model, and two thrusters and rudders are equivalent to one. Thus, the actuator model can be written as(7)τ=Bτactν,n,δ,where matrix B is the actuator configuration matrix and vector τactν,n,δ is the force and torques related to speed of the motor n and rudder angle δ:(8)Bν=20020−2lxR,τactν,n,δ=TL⊤,with(9)T=Tnnnn−Tnunu,L=Lδδ−Lδδδδuu.


## 2.2. Discrete Nonlinear Model

In the vessel experiments, the experiment data are collected in the discrete form. Also, vessel motion controllers are usually designed in the discrete zone. Therefore, the continuous-time model equation (1) is approximated as a discrete model, and the parameters of the discrete model are identified in this paper. The external force is all ignored for the convenience of discretizing the continuous-time model, and equation (1) can be rewritten as(10)Mν˙=Acα,where the matrix Ac and vector α denote the time-varying terms and constant parameters in the vessel model:(11)Acα=−Cνν−Dνν+Bτactν,n,δ,α=uuuu3r2vrn2nu000uvurvrvvvrrvrruuδuuδδuvurvrvvvrrvrruuδuuδδ⊤,Ac=XuXuuXuuumxg−12Nv˙+Yr˙m−Yv˙2Tnn−2Tnu0000Xu˙−mYvYr|Yvv|Yvr|YrvYrr2Lδ−2LδδYv˙−Xu˙12Nv˙+Yr˙−mxgNvNrNvvNvrNrvNrr−2lxRLδ2lxRLδδ⊤.The sample time of collecting data in simulation is denoted asTs; equation (10) can be approximated by the back Euler integration method at time k as(12)νk+1=νk+TsM−1Acαk=Bcαk+TsM−1Acαk,where Bc=100000000000100000000001000000⊤.Because the matrixBc and TsM−1Ac are both constant matrix, thus, equation (12) can be rewritten as(13)uk+1=a1uk+a2ukuk+a3uk3+a4rk2+a5vkrk+a6nk2+a7nkuk,(14)vk+1=b1ukvk+b2ukrk+b3vk+b4rk+b5vkvk+b6vkrk+b7rkvk+b8rkvk+b9ukukδk+b10ukukδkδk,(15)rk+1=c1ukvk+c2ukrk+c3vk+c4rk+c5vkvk+c6vkrk+c7rkvk+c8rkvk+c9ukukδk+c10ukukδkδk.The parametersa1,…,a7, b1,…,b10, and c1,…,c10 are unknown parameters. In Section 3, an ACO method is used to identify parameters of equations (13)–(15).


## 3. Parameter Identification Based on the ACO Method

In this section, the parameter of discrete vessel nonlinear model equations (13)–(15) is identified based on the ACO method.


### 3.1. Parameter Identification Formulation

Consider that the vessel has two main motion states: one is the straight-line surge motion and the other is the turning motion. Therefore, parameter identification includes two steps: the first step is to identify the parameters related to surge motion and the second step is to identify the residual parameters. In the straight-line surge motion, the lateral velocity and yaw angular velocity are both close to zero. Thus, the surge motion can be decoupled from the other motion, the parameters related to straight-line motion and thrust force can be identified. By ignoring the terms related to the lateral and yaw motion, equation (13) can be simplified as(16)uk+1=a1uk+a2ukuk+a6nk2+a7nkuk=θβuk,nk⊤,where θ=a1a2a6a7⊤ and βuk,nk⊤=ukukuknk2nkuk⊤.Based on equation (16), the predicted vessel surge velocity at time k+n, denoted as u^k+n+1, can be calculated by the motor speed at time k+n, denoted as nk+n, and predicted as surge velocity u^k+n. Then, the surge motion model identification problem can be described as a nonlinear optimization problem:(17)minfθ=12∑k=1Nu^k−uk2,s.t.u^1=u1,u^k=θβu^k−1,nk−1⊤,k≥2,where u^ denotes the predicted surge speed based on the parameter θ and u denotes the surge speed collected in the experiments.Since the parametersθ are identified in equation (17), then the rest parameters can also be described as a nonlinear optimization problem to identify. Based on equations (13)–(15), the vessel state can be predicted by(18)νk+1=θ2α1uk,vk,rk,nk,δk⊤+α2uk,nk⊤,where(19)θ2=000a4a500000b1b2b3b4b5b6b7b8b9b10c1c2c3c4c5c6c7c8c9c10,α2uk,nk⊤=θβuk,nk⊤00,α1=000r2vr00000uvurvrvvvrrvrruuδuuδδuvurvrvvvrrvrruuδuuδδ⊤.Based on equation (18), the predicted vessel state at time k+n, denoted as ν^k+n+1, can be calculated by the motor speed nk+n, rudder angle δk+n, the predicted surge velocity u^k+n, the predicted lateral velocity v^k+n, and the predicted yaw angular velocity r^k+n. Then, the parameters’ identification problem can be described as a nonlinear optimization problem:(20)minfθ2=12∑k=1Nν^k−νk,s.t.ν^1=ν1,ν^k+1=θ2α1u^k,v^k,r^k,nk,δk⊤+α2u^k,nk⊤,k≥2.The parameters now can be identified by solving two nonlinear optimization problem equations (17) and (20).


### 3.2. ACO Formulation

As mentioned above, the key to identify the model parameters is solving an optimization problem as equation (17) or equation (20). Generally, the typical unconstrained continuous optimization problem is shown in equation (21), where fx is the objective function and X is the decision variable (D.V) in the form of an N-dimensional vector in which the member changes continuously:(21)Minfx,X=x1x2⋯xN.As far as equations (17) and (20) are concerned, fθ and fθ2 are objective functions of the two identification problems, respectively. Accordingly, θ and θ2 are high-dimensional D.Vs.Traditional optimization algorithms find it difficult to solve the nonlinear high-dimensional optimization problems described in equations (17) and (20). A new computational paradigm called “Ant System” is proposed in [31], which is used for stochastic combinatorial optimization and motivated amounts of interrelated solutions to improve the performances. The representative one is Ant Colony Optimization (ACO). The parameter identification problem we focused on, however, is a typical optimization problem of which decision variables (D.V) vary in continuous domains.Despite the method that discretizes, the feasible region takes the effect to a certain extent; the accuracy and efficiency are severely restricted. Socha and Dorigo [32] extend ACO to continuous domains and keep the relevant conceptual structure invariable. The complete flowchart of the CACO is shown in Figure 1. In CACO, the fundamental work is initializing the elementary data structures and parameters to abstract the actual optimization problem and determine the initial state. Pheromone, inspired by the natural behavior of ants, which records the excellent solutions, is the significant feature of ACO. Different from the ACO which uses a high-dimensional matrix to represent the pheromone, CACO applies the continuous probability density function (CPDF), such as the Gaussian function. Furthermore, CACO maintains a solution archive of capacity k composed of kn-dimensional decision variables Sj1,Sj2,…,Sjn, which are constructed by ants, and the corresponding fitness fSj is to memorize the several former solutions. For the surge model identification problem, the variable n is 4 and the fitness is given by fθ. As for the remaining parameters identification problem, the variable n is determined by θ2, that is, 22, and the fitness is calculated by fθ2. Moreover, the solution vectors are supposed to be initialized using random values, sampling from the search space ℜ. The other related parameters which will be discussed below are supposed to be set reasonably at this stage.Figure 1
ACO method flowchart.As shown in Figure1, once the initialization work is accomplished, the algorithm will step into searching for solutions in ℜ and further instruct to update the pheromone. In the course of searching for solutions for an individual ant, the first step is choosing a target; that is, select a solution from the archive based on the selected probability following equations (22) and (23) and the roulette method. The value l is the rank of the corresponding solution, k is the size of the archive, and q is an extra parameter to adjust the standard deviation, which has an effect on the shot probability of different ranks. Namely, q is a trade-off factor between the global and local optimal solutions:(22)ω1=1qk2πe−l−12/2q2k2,(23)pl=ω1∑r=1kωr.It is worth noting that the ants in ACO or CACO are not a container (i.e., a vector composed of decision variables) but a solution builder (i.e., a pathfinder in the search space to construct the solution vector). For ann-dimensional decision variable, an ant needs n steps to accomplish the construction per iteration. Once the target is selected, the destination of the ith dimension variable is assigned by Gaussian sampling using the function described as(24)glix,μli,σli=1σli2πe−x−μli/2μli2.The parameterμli=sli, and σli is calculated as equation (25), where ξ has equivalent effect with the pheromone evaporation rate which can influence the speed of convergence:(25)σli=ξ∑j=1ksji−slik−1.On account of the pheromone reflects on the solution archive, an update-pheromone operation may be accomplished by renewing the table. After all of the ants complete the construction of solutions, there will bem+k solutions. For the sake of the superior solution to maintain the archive, all solutions should be ordered according to the evaluation index (e.g., the root mean square error between predicting results and sample data). In the subsequent step, the worst m solutions would be discarded and the others remain in the archive.After updating the archive, the algorithm tests’ termination conditions (such as the maximum iterations and convergence precision) are in accordance with the optimization results of the current generation. If the termination condition is not met, the optimization problem will repeat the searching, evaluating, and renewing process mentioned above until the condition is satisfied. Otherwise, the iteration should be broken up and the optimization results should be output, which are the optimized model parametersθ and θ2 and the corresponding model performance.


## 3.1. Parameter Identification Formulation

Consider that the vessel has two main motion states: one is the straight-line surge motion and the other is the turning motion. Therefore, parameter identification includes two steps: the first step is to identify the parameters related to surge motion and the second step is to identify the residual parameters. In the straight-line surge motion, the lateral velocity and yaw angular velocity are both close to zero. Thus, the surge motion can be decoupled from the other motion, the parameters related to straight-line motion and thrust force can be identified. By ignoring the terms related to the lateral and yaw motion, equation (13) can be simplified as(16)uk+1=a1uk+a2ukuk+a6nk2+a7nkuk=θβuk,nk⊤,where θ=a1a2a6a7⊤ and βuk,nk⊤=ukukuknk2nkuk⊤.Based on equation (16), the predicted vessel surge velocity at time k+n, denoted as u^k+n+1, can be calculated by the motor speed at time k+n, denoted as nk+n, and predicted as surge velocity u^k+n. Then, the surge motion model identification problem can be described as a nonlinear optimization problem:(17)minfθ=12∑k=1Nu^k−uk2,s.t.u^1=u1,u^k=θβu^k−1,nk−1⊤,k≥2,where u^ denotes the predicted surge speed based on the parameter θ and u denotes the surge speed collected in the experiments.Since the parametersθ are identified in equation (17), then the rest parameters can also be described as a nonlinear optimization problem to identify. Based on equations (13)–(15), the vessel state can be predicted by(18)νk+1=θ2α1uk,vk,rk,nk,δk⊤+α2uk,nk⊤,where(19)θ2=000a4a500000b1b2b3b4b5b6b7b8b9b10c1c2c3c4c5c6c7c8c9c10,α2uk,nk⊤=θβuk,nk⊤00,α1=000r2vr00000uvurvrvvvrrvrruuδuuδδuvurvrvvvrrvrruuδuuδδ⊤.Based on equation (18), the predicted vessel state at time k+n, denoted as ν^k+n+1, can be calculated by the motor speed nk+n, rudder angle δk+n, the predicted surge velocity u^k+n, the predicted lateral velocity v^k+n, and the predicted yaw angular velocity r^k+n. Then, the parameters’ identification problem can be described as a nonlinear optimization problem:(20)minfθ2=12∑k=1Nν^k−νk,s.t.ν^1=ν1,ν^k+1=θ2α1u^k,v^k,r^k,nk,δk⊤+α2u^k,nk⊤,k≥2.The parameters now can be identified by solving two nonlinear optimization problem equations (17) and (20).


## 3.2. ACO Formulation

As mentioned above, the key to identify the model parameters is solving an optimization problem as equation (17) or equation (20). Generally, the typical unconstrained continuous optimization problem is shown in equation (21), where fx is the objective function and X is the decision variable (D.V) in the form of an N-dimensional vector in which the member changes continuously:(21)Minfx,X=x1x2⋯xN.As far as equations (17) and (20) are concerned, fθ and fθ2 are objective functions of the two identification problems, respectively. Accordingly, θ and θ2 are high-dimensional D.Vs.Traditional optimization algorithms find it difficult to solve the nonlinear high-dimensional optimization problems described in equations (17) and (20). A new computational paradigm called “Ant System” is proposed in [31], which is used for stochastic combinatorial optimization and motivated amounts of interrelated solutions to improve the performances. The representative one is Ant Colony Optimization (ACO). The parameter identification problem we focused on, however, is a typical optimization problem of which decision variables (D.V) vary in continuous domains.Despite the method that discretizes, the feasible region takes the effect to a certain extent; the accuracy and efficiency are severely restricted. Socha and Dorigo [32] extend ACO to continuous domains and keep the relevant conceptual structure invariable. The complete flowchart of the CACO is shown in Figure 1. In CACO, the fundamental work is initializing the elementary data structures and parameters to abstract the actual optimization problem and determine the initial state. Pheromone, inspired by the natural behavior of ants, which records the excellent solutions, is the significant feature of ACO. Different from the ACO which uses a high-dimensional matrix to represent the pheromone, CACO applies the continuous probability density function (CPDF), such as the Gaussian function. Furthermore, CACO maintains a solution archive of capacity k composed of kn-dimensional decision variables Sj1,Sj2,…,Sjn, which are constructed by ants, and the corresponding fitness fSj is to memorize the several former solutions. For the surge model identification problem, the variable n is 4 and the fitness is given by fθ. As for the remaining parameters identification problem, the variable n is determined by θ2, that is, 22, and the fitness is calculated by fθ2. Moreover, the solution vectors are supposed to be initialized using random values, sampling from the search space ℜ. The other related parameters which will be discussed below are supposed to be set reasonably at this stage.Figure 1
ACO method flowchart.As shown in Figure1, once the initialization work is accomplished, the algorithm will step into searching for solutions in ℜ and further instruct to update the pheromone. In the course of searching for solutions for an individual ant, the first step is choosing a target; that is, select a solution from the archive based on the selected probability following equations (22) and (23) and the roulette method. The value l is the rank of the corresponding solution, k is the size of the archive, and q is an extra parameter to adjust the standard deviation, which has an effect on the shot probability of different ranks. Namely, q is a trade-off factor between the global and local optimal solutions:(22)ω1=1qk2πe−l−12/2q2k2,(23)pl=ω1∑r=1kωr.It is worth noting that the ants in ACO or CACO are not a container (i.e., a vector composed of decision variables) but a solution builder (i.e., a pathfinder in the search space to construct the solution vector). For ann-dimensional decision variable, an ant needs n steps to accomplish the construction per iteration. Once the target is selected, the destination of the ith dimension variable is assigned by Gaussian sampling using the function described as(24)glix,μli,σli=1σli2πe−x−μli/2μli2.The parameterμli=sli, and σli is calculated as equation (25), where ξ has equivalent effect with the pheromone evaporation rate which can influence the speed of convergence:(25)σli=ξ∑j=1ksji−slik−1.On account of the pheromone reflects on the solution archive, an update-pheromone operation may be accomplished by renewing the table. After all of the ants complete the construction of solutions, there will bem+k solutions. For the sake of the superior solution to maintain the archive, all solutions should be ordered according to the evaluation index (e.g., the root mean square error between predicting results and sample data). In the subsequent step, the worst m solutions would be discarded and the others remain in the archive.After updating the archive, the algorithm tests’ termination conditions (such as the maximum iterations and convergence precision) are in accordance with the optimization results of the current generation. If the termination condition is not met, the optimization problem will repeat the searching, evaluating, and renewing process mentioned above until the condition is satisfied. Otherwise, the iteration should be broken up and the optimization results should be output, which are the optimized model parametersθ and θ2 and the corresponding model performance.


## 4. Simulation Results

In this paper, the CyberShip II is used as a study example. The kinetics model structure of the vessel is already known, and the parameters of the kinetics model are unknown. The parameters are identified by the ACO method and the LS method, respectively, and the identification results are compared to demonstrate the effectiveness of the ACO method. Normally, in the vessel identification experiments, the vessel states and actuator states such as velocity, angular velocity, motor speed, and rudder angle are collected by the sensor with noise. Therefore, the Gaussian white noise is added in the whole identification simulations to simulate the sensor noise. The whole identification procedure is divided into two steps, the first step is to identify the parameters related to the surge motion and the second step is to identify the parameters related to the sway and yaw motion.


### 4.1. Surge Model Identification

In the vessel surge motion, sway velocity and yaw angular velocity are approximately equal to zero. Therefore, the couple terms of the vessel kinetics model related to sway and yaw motion can be ignored. The purpose of the surge model identification is to identify the parametersa1, a2, a3, a5, and a6 in equation (13).The vessel motor speed is used as the excitation input which varies randomly between 5 rpm and 25 rpm, and the vessel surge velocity as the response to the excitation input varies between 0.2m/s and 0.8 m/s, as shown in Figure 2. The transient motor speed and surge speed are used as training data to identify the parameters. The parameters related to the ACO method in solving surge model identification problem are shown Table 1, and the ACO method convergence curve are shown in Figure 3. From Figure 3, we can see the fitness which is the value of fθ in equation (17) gradually decreases and becomes stable after 450 iterations. The identified results of the ACO method and LS method are shown in Table 2.Figure 2
Identification data of the surge model.
(a)(b)Table 1
Parameters of the ACO method in surge model identification.
ParametersnmkqξGValue520101e − 030.75500Figure 3
ACO-method convergence curve.Table 2
Parameters of the surge model.
Parametersa1a2a3a6a7LS0.9953−0.0060−0.03882.6357e − 055.6225e − 04ACO0.9972−0.0051−0.02272.8218e − 05−2.8661e − 06In Figure2, the prediction results of parameters identified by the different methods are very close. In addition, another two different surge motions are used as validation data to compare the performance of different algorithms, as shown in Figures 4 and 5. In Figure 4, the input motor speed is a PRBS sequence. From the results, we can see the prediction results of different methods are very close, except that the prediction results of the LS method are slower than the simulation data when the surge speed is larger than 0.9 m/s.Figure 4
Validation data 1 (results of the surge model under PRBS motor speed).
(a)(b)Figure 5
Validation data 2 (results of the surge model under different constant motor speeds).
(a)(b)In Figure5, the input motor speed is maintained at 5 rpm, 10 rpm, and 20 rpm, so the surge speed can finally be maintained at around 0.16 m/s, 0.36 m/s, and 0.67 m/s. We can see from the results that the prediction results of different methods are almost the same when the motor speed is kept at 5 rpm and 10 rpm. However, the prediction surge speed of the LS method is larger than the simulation results when the motor speed is kept at 10 rpm.In order to quantify the comparison result of different methods, theR-square is used as an index function. The closer the R-square to 1, the better the identified results, and in some situations, the R-square may be negative; then, we use the symbol − to indicate it. From the R-square results of different methods in Table 3, we can see the R-square values of the ACO method and the LS method are very close. However, the R-square value of the ACO method is better than the LS method in training data and validation data, so the ACO method has better performance than the LS method in the surge model parameters’ identification simulation experiments.Table 3
R-square values of different methods on the surge model.
R-squareTraining dataValidation data 1Validation data 2LS method0.99530.99590.9908ACO method0.99590.99820.9934


### 4.2. Lateral and Yaw Model Identification

In the sway and yaw model identification simulation experiments, the rest parameters related to the sway and yaw motion in equations (13)–(15) are identified. Zigzag tests are widely used in vessel parameters’ identification experiments [16, 24]. In this paper, a 20°/20° zigzag test contains constant motor speed which is used to generate the training data, as shown in Figure 6 and Table 4.Figure 6
Identification results of vessel dynamic.
(a)(b)(c)(d)(e)(f)Table 4
Parameters of the ACO method in lateral and yaw model identification.
ParametersnmkqξGValue2230151e− 030.751000The convergence process of the ant colony algorithm is shown in Figure7. The result shows that the optimal solution almost converges to the expected range after 1000 generations. The parameters identified by the LS method and ACO method are shown in Table 5.Figure 7
ACO method convergence curve.Table 5
Parameters of lateral and yaw motion.
Parametersa3a4b1b2LS0.00240.1292−3.6851e− 04−0.0789ACO0.00410.13020.0012−0.0780Parametersb3b4b5b6Ls0.9978−0.0206−0.1089−3.1260e− 4ACO0.9917−0.0262−0.1025−0.0030Parametersb7b8b9b10Ls−0.0026−0.01260.02380.0083ACO−0.0019−0.00900.0443−0.0287Parametersc1c2c3c4Ls−0.04130.0013−1.3609e− 040.9905ACO−0.0453−0.00133.2504e− 040.9864Parametersc5c6c7c8Ls0.0277−0.00170.0011−8.6515e− 04ACO0.02035.3517e− 047.6613e− 04−0.0039Parametersc9c10Ls−0.0243−0.0085ACO−0.03420.0251From Figure6 we can see that the LS method has better prediction performance than the ACO method at the surge model, sway model, and yaw model. To further verify the effectiveness and the performance of the ACO method, we use another two different excitation inputs to generate the validation data which are shown in Figures 8 and 9, and for the convenience of comparison, the heading angle is not changed to the range from 0° to 360°. In addition, the R-square value is used to quantitatively compare the prediction performance for the two different methods, as shown Table 6.Figure 8
Validation results of vessel dynamic by the least-square method.
(a)(b)(c)(d)(e)(f)Figure 9
Validation results of vessel dynamic by the least-square method.
(a)(b)(c)(d)(e)(f)Table 6
R-square values of different methods on the vessel model.
Identification errorSurge modelSway modelYaw modelHeading modelLS0.98900.99570.99520.9964ACO0.84780.91820.94170.6804LS0.91060.94580.97400.9955ACO0.94760.96720.98310.9982LS0.7050———ACO0.95260.97450.97840.9942In Figure8, the motor speed maintains at 10 rpm, 20 rpm, and 30 rpm at different times, and the rudder angle is a PRBS sequence with different amplitudes. From the results, we can see the prediction results of ACO are slightly better than the LS method at a motor speed of 10 rpm, but the prediction results of the two methods are similar at a motor speed of 20 rpm and 30 rpm. The R-square value of the ACO method is closer to 1 than the LS method on the surge model, sway model, yaw model, and heading model.In Figure9, the motor speed maintains at 10 rpm, 20 rpm, and 30 rpm at different times, and the rudder angle maintains at 5°, −15°, and 25° at different times. From the results, we can see only at a motor speed of 5 rpm and rudder angle of 5°; the two methods have similar prediction results. From the R-square values, we can see the prediction results of the ACO method have very good performance, but the prediction results of the LS method are only valid at low motor speed and small rudder angle.In the identification results of the sway and yaw motion, although the prediction results of the LS method are better than the LS method in training data, the identified parameters of the LS method are invalid in predicting the surge, sway, and yaw motion when the motor speed and rudder angle maintain the constant value. However, the parameters identified by the ACO method have a good prediction result in the surge, sway, and yaw motion both in the training data and validation data.


## 4.1. Surge Model Identification

In the vessel surge motion, sway velocity and yaw angular velocity are approximately equal to zero. Therefore, the couple terms of the vessel kinetics model related to sway and yaw motion can be ignored. The purpose of the surge model identification is to identify the parametersa1, a2, a3, a5, and a6 in equation (13).The vessel motor speed is used as the excitation input which varies randomly between 5 rpm and 25 rpm, and the vessel surge velocity as the response to the excitation input varies between 0.2m/s and 0.8 m/s, as shown in Figure 2. The transient motor speed and surge speed are used as training data to identify the parameters. The parameters related to the ACO method in solving surge model identification problem are shown Table 1, and the ACO method convergence curve are shown in Figure 3. From Figure 3, we can see the fitness which is the value of fθ in equation (17) gradually decreases and becomes stable after 450 iterations. The identified results of the ACO method and LS method are shown in Table 2.Figure 2
Identification data of the surge model.
(a)(b)Table 1
Parameters of the ACO method in surge model identification.
ParametersnmkqξGValue520101e − 030.75500Figure 3
ACO-method convergence curve.Table 2
Parameters of the surge model.
Parametersa1a2a3a6a7LS0.9953−0.0060−0.03882.6357e − 055.6225e − 04ACO0.9972−0.0051−0.02272.8218e − 05−2.8661e − 06In Figure2, the prediction results of parameters identified by the different methods are very close. In addition, another two different surge motions are used as validation data to compare the performance of different algorithms, as shown in Figures 4 and 5. In Figure 4, the input motor speed is a PRBS sequence. From the results, we can see the prediction results of different methods are very close, except that the prediction results of the LS method are slower than the simulation data when the surge speed is larger than 0.9 m/s.Figure 4
Validation data 1 (results of the surge model under PRBS motor speed).
(a)(b)Figure 5
Validation data 2 (results of the surge model under different constant motor speeds).
(a)(b)In Figure5, the input motor speed is maintained at 5 rpm, 10 rpm, and 20 rpm, so the surge speed can finally be maintained at around 0.16 m/s, 0.36 m/s, and 0.67 m/s. We can see from the results that the prediction results of different methods are almost the same when the motor speed is kept at 5 rpm and 10 rpm. However, the prediction surge speed of the LS method is larger than the simulation results when the motor speed is kept at 10 rpm.In order to quantify the comparison result of different methods, theR-square is used as an index function. The closer the R-square to 1, the better the identified results, and in some situations, the R-square may be negative; then, we use the symbol − to indicate it. From the R-square results of different methods in Table 3, we can see the R-square values of the ACO method and the LS method are very close. However, the R-square value of the ACO method is better than the LS method in training data and validation data, so the ACO method has better performance than the LS method in the surge model parameters’ identification simulation experiments.Table 3
R-square values of different methods on the surge model.
R-squareTraining dataValidation data 1Validation data 2LS method0.99530.99590.9908ACO method0.99590.99820.9934


## 4.2. Lateral and Yaw Model Identification

In the sway and yaw model identification simulation experiments, the rest parameters related to the sway and yaw motion in equations (13)–(15) are identified. Zigzag tests are widely used in vessel parameters’ identification experiments [16, 24]. In this paper, a 20°/20° zigzag test contains constant motor speed which is used to generate the training data, as shown in Figure 6 and Table 4.Figure 6
Identification results of vessel dynamic.
(a)(b)(c)(d)(e)(f)Table 4
Parameters of the ACO method in lateral and yaw model identification.
ParametersnmkqξGValue2230151e− 030.751000The convergence process of the ant colony algorithm is shown in Figure7. The result shows that the optimal solution almost converges to the expected range after 1000 generations. The parameters identified by the LS method and ACO method are shown in Table 5.Figure 7
ACO method convergence curve.Table 5
Parameters of lateral and yaw motion.
Parametersa3a4b1b2LS0.00240.1292−3.6851e− 04−0.0789ACO0.00410.13020.0012−0.0780Parametersb3b4b5b6Ls0.9978−0.0206−0.1089−3.1260e− 4ACO0.9917−0.0262−0.1025−0.0030Parametersb7b8b9b10Ls−0.0026−0.01260.02380.0083ACO−0.0019−0.00900.0443−0.0287Parametersc1c2c3c4Ls−0.04130.0013−1.3609e− 040.9905ACO−0.0453−0.00133.2504e− 040.9864Parametersc5c6c7c8Ls0.0277−0.00170.0011−8.6515e− 04ACO0.02035.3517e− 047.6613e− 04−0.0039Parametersc9c10Ls−0.0243−0.0085ACO−0.03420.0251From Figure6 we can see that the LS method has better prediction performance than the ACO method at the surge model, sway model, and yaw model. To further verify the effectiveness and the performance of the ACO method, we use another two different excitation inputs to generate the validation data which are shown in Figures 8 and 9, and for the convenience of comparison, the heading angle is not changed to the range from 0° to 360°. In addition, the R-square value is used to quantitatively compare the prediction performance for the two different methods, as shown Table 6.Figure 8
Validation results of vessel dynamic by the least-square method.
(a)(b)(c)(d)(e)(f)Figure 9
Validation results of vessel dynamic by the least-square method.
(a)(b)(c)(d)(e)(f)Table 6
R-square values of different methods on the vessel model.
Identification errorSurge modelSway modelYaw modelHeading modelLS0.98900.99570.99520.9964ACO0.84780.91820.94170.6804LS0.91060.94580.97400.9955ACO0.94760.96720.98310.9982LS0.7050———ACO0.95260.97450.97840.9942In Figure8, the motor speed maintains at 10 rpm, 20 rpm, and 30 rpm at different times, and the rudder angle is a PRBS sequence with different amplitudes. From the results, we can see the prediction results of ACO are slightly better than the LS method at a motor speed of 10 rpm, but the prediction results of the two methods are similar at a motor speed of 20 rpm and 30 rpm. The R-square value of the ACO method is closer to 1 than the LS method on the surge model, sway model, yaw model, and heading model.In Figure9, the motor speed maintains at 10 rpm, 20 rpm, and 30 rpm at different times, and the rudder angle maintains at 5°, −15°, and 25° at different times. From the results, we can see only at a motor speed of 5 rpm and rudder angle of 5°; the two methods have similar prediction results. From the R-square values, we can see the prediction results of the ACO method have very good performance, but the prediction results of the LS method are only valid at low motor speed and small rudder angle.In the identification results of the sway and yaw motion, although the prediction results of the LS method are better than the LS method in training data, the identified parameters of the LS method are invalid in predicting the surge, sway, and yaw motion when the motor speed and rudder angle maintain the constant value. However, the parameters identified by the ACO method have a good prediction result in the surge, sway, and yaw motion both in the training data and validation data.


## 5. Conclusions

This paper presents the parameters’ identification under the known structure of the vessel kinetics model. The parameters’ problem is summarized as a least-square problem; the least-square method and the ant colony optimization method are used to solve the least-square problem. The solution to the least-square problem is the parameters to be identified. The identification procedure is divided into two parts. The first part of the identification procedure is to identify the parameters related to surge motion. The second part of the identification procedure is to identify the rest parameters of the vessel kinetics model. In the identification procedure, the transient excitation inputs are used to generate the training data; the constant excitation inputs and different transient inputs are used to generate the validation data. TheR-square value is used as an index function to quantitatively compare the prediction results of the ACO method and the LS method. The comparison of the identification results and maneuvering predictions demonstrate the effectiveness of the ACO method and reflect the performance advantage of the ACO method.



---

*Source: 6256785-2021-07-26.xml*