2. Nonlinear ship model

Ship motion can be described by a state space model or an input-output model (Zhang et al., Reference Zhang and Zhang2020a, Reference Zhang, Zhang and Ren2020b, Reference Zhang2020c). The former can deal with the multiple variables associated with ship motion. Thus, introducing wind, waves and current disturbances can generate more accurate and direct results. However, it comes at the price of greater computational complexity. The latter approach is also known as a response model (Deng and Zhang, Reference Deng and Zhang2020). It ignores drift speed and focuses on the main line of the dynamic rudder angle, the yaw angular velocity and the heading angle $(\delta \to \dot{\psi } \to \psi )$. The resulting differential equations can retain factors of nonlinear influence. This approach can even convert the disturbances arising from the wind and waves into a disturbance rudder angle, which can then form an input signal together with the actual rudder angle (see Figure 1). This model is a generalisation of the linear Nomoto model (Nomoto et al., Reference Nomoto, Taguchi, Honda and Hirnao1957; Zhang and Jin, Reference Zhang, Yao, Xu and Zhang2013).

Figure 1. Nonlinear ship model

This paper adopts a response ship mathematical model as shown in Figure 1, which is composed of a first-order Nomoto model and a nonlinear feedback compensating item. According to marine practice, a set of rudder servo system is considered, angle limiter and revolution rate limiter. Rudder angle $\delta \in [{ - 35^\circ ,35^\circ } ]$. The first-order Nomoto model from $\delta$ to yaw rate r is presented as

(1)\begin{equation}{G_{r\delta }}(s) = \frac{K}{{Ts + 1}}\end{equation}

where the nonlinear feedback compensating item $f(u)$ is expressed as

(2)\begin{equation}f(u) = (\alpha - 1/K)\dot{\psi } + \beta {\dot{\psi }^3}\end{equation}

The nonlinear Nomoto model is then used as the ship mathematical model (Zhang and Jin, Reference Zhang, Yao, Xu and Zhang2013). The relevant equation is:

(3)\begin{equation}\ddot{\psi } + \frac{K}{T}(\alpha \mathop {\dot{\psi } + \beta {{\dot{\psi }}^3}}\limits^{} )=\frac{K}{T}\delta\end{equation}

where $\psi $ is the heading angle; $\dot{\psi }$ is the yaw rate; $\ddot{\psi }$ is the angular acceleration;$\delta $ is the rudder angle; K and T are ship manoeuvrability indices; and $\alpha $ and $\beta $ are the nonlinear coefficients of the yaw rate $\dot{\psi }$.

When the ship makes a steady turning test, the mathematical model shown in Equation (3) can be further simplified. If we let $r = \dot{\psi }$, $\dot{r} = \ddot{\psi }$ and $\dot{r} = \ddot{\psi } = 0$ after the turning is stable (Jia and Zhang, Reference Jia and Zhang2001), the equation can be rewritten as:

(4)\begin{equation}\delta =\alpha \mathop {\dot{\psi } + \beta {{\dot{\psi }}^3}}\limits^{} =\alpha r + \beta {r^3}=\left[ {\begin{array}{*{20}{c}} r&{{r^3}} \end{array}} \right]{\left[ {\begin{array}{*{20}{c}} \alpha &\beta \end{array}} \right]^{\rm T}}\end{equation}

Note, here, that the input and output of Equation (4) and Equation (3) are opposite. For Equation (3), $\delta $ is the input, but it is the output for Equation (4). In this paper, Equation (4) is used as the mathematical ship model so as to be able to identify the nonlinear coefficients $\alpha ,\beta $.

3. Nonlinear innovation stochastic gradient identification algorithm

The linear regression model is:

(5)\begin{equation}y(t) = {\varphi^{\rm T}}\theta + v(t)\end{equation}

where $y(t)$ is the output. In Equation (2), this is $\delta $. $\varphi (t)$ is the regression information vector. In Equation (4), it is ${\left[ {\begin{array}{*{20}{c}} r&{{r^3}} \end{array}} \right]^{\rm T}}$. $\theta $ relates to the parameters to be identified. In Equation (5), they are ${\left[ {\begin{array}{*{20}{c}} \alpha &\beta \end{array}} \right]^{\rm T}}$.$v(t)$ is zero-mean-value random noise. For the identification system described by Equation (5), the stochastic gradient (SG) will be

(6)\begin{align}\hat{\theta }(t) &= \hat{\theta }(t - 1) + \frac{{\varphi (t)}}{{r(t)}}e(t)\end{align}

(7)\begin{align}e(t) &= y(t) - {\varphi ^{\rm T}}(t)\hat{\theta }(t - 1)\end{align}

(8)\begin{align}r(t) &= r(t - 1) + {||{\varphi (t)} ||^2}\end{align}

where $\hat{\theta }(t)$ and $\hat{\theta }(t - 1)$ are estimates for $\theta $ at the current and previous moment, respectively. $e(t)$ is innovation, which means useful information that can improve the accuracy of the parameter estimation or state estimation. Unlike the least squares method, the SG algorithm does not need to compute a covariance matrix and has less computational complexity. However, its convergence speed is slower and its identification is less accurate. Drawing upon the nonlinear feedback algorithm of Zhang et al. (2016, Reference Zhang and Jin2018) and the multi-innovation identification algorithm of Ding (Reference Ding2013), the innovation can be processed by a nonlinear arc tangent function. The improved SG algorithm based on nonlinear innovation (NISG) can therefore be written as follows:

(9)\begin{align}\hat{\theta }(t) &= \hat{\theta }(t - 1) + \frac{{\varphi (t)}}{{r(t)}}{{\rm e}^{\prime}}(t)\end{align}

(10)\begin{align}{{\rm e}^{\prime}}(t) &= A\arctan ({\omega _{}}e(t)) = \arctan (\omega (y(t) - {\varphi ^{\rm T}}(t)\hat{\theta }(t - 1)))\end{align}

(11)\begin{align}r(t) &= r(t - 1) + {||{\varphi (t)} ||^2}\end{align}

In Equation (10), the value of A, $\omega $ is a nonlinear parameter, generally 0.5–2.0.

The performance of this algorithm can be analysed using the Martingale convergence theorem and stochastic process theory (Ding and Yang, Reference Ding and Yang1999). Let $\{{v(t)} \}$ be a Martingale difference sequence defined in the probability space {K,F,P}, where {F_t} is an algebraic sequence generated by $\{{v(t)} \}$, and $\{{v(t)} \}$ satisfies the following noise hypothesis:

(12)\begin{equation}\begin{aligned} &1)\;E[v(t)|{F_{t - 1}}] = 0,a.s.;\\ &2)\;E[{v^2}(t)|{F_{t - 1}}] = e_v^2(t) \le e_v^2 \lt \infty ,a.s.;\\ &3)\;\mathop {\lim \sup }\limits_{t \to \infty } \dfrac{1}{t}\sum\limits_{t = 1}^t {{v^2}(i)} \le e_v^2 \lt \infty ,a.s.. \end{aligned}\end{equation}

Note that a.s. here means ‘almost surely’.

Lemma 1: For the SG algorithms (9) and (11), indicates that the following conclusions are valid:

(13)\begin{equation}\begin{aligned} &{\rm a})\;\mathop {\lim }\limits_{t \to \infty } \sum\limits_{i = 0}^t {\dfrac{{||\varphi (i)|{|^2}}}{{r(i - 1)r(i)}}} \lt \infty ,a.s.;\\ &{\rm b)}\;e(t) = \dfrac{{r(t - 1)}}{{r(t)}}Z(t). \end{aligned}\end{equation}

where

\[\begin{split} e(t) &= y(t) - {\varphi ^{\rm T}}(t)\hat{\theta }(t)\\ Z(t) &= y(t) - {\varphi ^{\rm T}}(t)\hat{\theta }(t - 1) \end{split}\]

Lemma 2: For the same algorithms states that, if there are constants T ₀, U and the integer N, that can establish the following strong persistent excitation conditions:

(14)\begin{equation}\begin{split} &4)\;{T_0}I \le \dfrac{1}{N}\sum\limits_{I = 0}^{N - 1} {\varphi (t - i){\varphi ^{\rm T}}} (t - i) \le UI,a.s.\\ &\quad 0 \lt {T_0} \le U \lt \infty ,N \ge n. \end{split}\end{equation}

then, $r(t)$ satisfies:

\[\begin{split} r(t - N) + nN{T_0} \le r(t) \le r(t - N) + nNU;\\ n{T_0}(t - N) \le r(t) \le nU(t + N) + 1;\\ 0 \lt n{T_0} \le \mathop {\lim }\limits_{t \to \infty } \dfrac{{r(t)}}{t} \le nU \lt \infty ,a.s.. \end{split}\]

If the noise hypotheses 1)–3) in Equation (12) and the strong persistent excitation condition 4) in Equation (14) hold, then the estimated parameters given by the SG algorithms Equations (6) and (8) almost surely converge to the true parameters, i.e., $\mathop {\lim \hat{\theta }}\limits_{t \to \infty } (t) = \theta ,a.s.$. The specific proof process for this can be found in Ding and Yang (Reference Ding and Yang1999).

For the improved SG algorithm based on nonlinear innovation, $\mathrm{e^{\prime}}(t) = \arctan ({\omega _{}}e(t))$ in Equation (8) can be expressed by a Taylor series expansion:

(15)\begin{equation}\arctan (\omega e(t)) \approx \omega e(t) - \frac{{{{(\omega e(t))}^3}}}{3}{\rm + }\frac{{{{(\omega e(t))}^5}}}{5} - \frac{{{{(\omega e(t))}^7}}}{7}\ldots \end{equation}

$\omega e(t)$ in Equation (15) is generally small in marine practice. So, $\arctan (\omega e(t)) \approx \omega e(t)$. In the above equations, $r(t)$ then becomes $r^{\prime}(t)=A\omega r(t)$ and Equations (12)–(14) still hold. Thus, the estimated parameters given by the improved SG algorithms (9) and (10) almost surely converge to the true parameters.

4. Simulation and results

This paper uses Dalian Maritime University's new teaching and training ship Yu Peng to test an application of the ship mathematical model. Yu Peng made its first voyage in 2017, so its test data are relatively complete. This helps with the verification of the identification of the ship model parameters. The ship mathematical model is based on the Nomoto model presented in Equation (1). The ship parameters required for establishing the Nomoto model are shown in Table 1. The experimental simulation was carried out using Visual Basic.

Table 1. Particulars of Yu Peng

According to the particulars of Yu Peng given in Table 1, the true values of the Nomoto model parameters for Yu Peng were calculated in Visual Basic, as shown in Table 2 (Nomoto et al., Reference Nomoto, Taguchi, Honda and Hirnao1957; Zhang and Zhang, Reference Zhang, Zhang and Pang2019). On the basis of Table 2, a Visual Basic program was designed to carry out experimental simulations and produce identification data. In the simulations, random interference with an amplitude of 0.1 for $\delta $ and random interference with an amplitude of 0.5 for ψ were introduced, according to marine practice. Twenty-six sets of turning test data were then generated, from hard port $({ - 35^\circ } )$ to hard starboard $({35^\circ } )$ in steps of $2.5^\circ $.

Table 2. Mathematical model parameters for Yu Peng

Next, the identification capacity of the least squares (LS) method, SG algorithm, and the NISG algorithm were compared for the nonlinear parameters, $\alpha $,$\beta $, across just the 26 sets of test data. The identification error for parameters $\alpha $,$\beta $ means the ratio of the difference (identification value and true value) and the true value.

As shown in Figure 2, the nonlinear parameters $\alpha $,$\beta $ were identified by LS and SG. The identification result for LS is the solid line and indicates an error of 14.01% for $\alpha $ and an error of 19.74% for $\beta $, with the mean error being 16.86%. The identification result for SG (the dashed line) indicates an error of 91.54% for $\alpha $ and an error of over 100% for $\beta $. So, on the basis of only 26 sets of test data, the error for SG is too large to be suitable for parameter identification, which means it has no practical significance. LS still maintains some level of identification accuracy when there is less test data and it is obviously better than SG.

Figure 2. Comparison of the identification results (LS and SG) for Yu Peng

For the improved SG algorithm based on nonlinear innovation, the refined models for identifying $\alpha $ and $\beta $ are as follows:

(16)\begin{align}\hat{\alpha }(t) &= \hat{\alpha }(t - 1) + \frac{{\varphi (t)}}{{r(t)}}\arctan (0.5\ast ((y(t) - {\varphi ^{\rm T}}(t)\hat{\alpha }(t - 1)))\end{align}

(17)\begin{align}\hat{\beta }(t) &= \hat{\beta }(t - 1) + \frac{{\varphi (t)}}{{r(t)}}\arctan 1.2\mathrm{\ast }(0.9\ast ((y(t) - {\varphi ^{\rm T}}(t)\hat{\beta }(t - 1)))\end{align}

As can be seen in Figure 3, the nonlinear parameters $\alpha $ and $\beta $ were also identified by NISG. The identification result for LS (the solid line) is reproduced for easier comparison. As previously stated, the error for $\alpha $ was 14.01% and the error for $\beta $ was 19.74%, with a mean error of 16.86%. The identification result for NISG (the dashed line) indicates an error of 0.23% for $\alpha $ and an error of 9.19% for $\beta $, with a mean error of only 4.71%. So, the identification accuracy for NISG can reach 95.3%.

Figure 3. Comparison of the identification results (LS and NISG) for Yu Peng

Based on just 26 sets of test data, the results suggest that the NISG method has a 12.2% higher identification accuracy than the LS method. It can also be seen in Figure 4 that the curve is smooth, with less fluctuation, and the identification speed is faster than LS.

Figure 4. Comparison of identification results (LS and NISG) for Yu Kun

To further verify the effectiveness of the NISG method, another of Dalian Maritime University's teaching and training ships, the Yu Kun, was also used as the basis of an experimental simulation. The parameters of the Yu Kun are shown in Tables 3 and 4.

Table 3. Particulars of Yu Kun

Table 4. Mathematical model parameters for Yu Kun

Again, the nonlinear parameters $\alpha $ and $\beta $ of the Yu Kun were identified by LS and NISG. The models described by Equations (14) and (15) were used once more for NISG.

As shown in Figure 4, the identification results for LS (the solid line) indicate an error of 12.08% for $\alpha $ and an error for $\beta $ of 23.27%, with a mean error of 17.68%. The results for NISG (the dashed line) show that the error for $\alpha $ was 9.44% and the error for $\beta $ was 3.08%, with a mean error of only 6.26%. So, it achieved an identification accuracy of about 93.7%.

For the ship Yu Kun, the identification accuracy of the NISG method was still 11.4% higher than the LS method under the circumstances of only a small amount of test data, thus proving the validity and generalisability of the algorithm.