Control for Ship Course-Keeping Using Optimized Support Vector Machines

Abstract

Support vector machines (SVM) are proposed in order to obtain a robust controller for ship course-keeping. A cascaded system is constructed by combining the dynamics of the rudder actuator with the dynamics of ship motion. Modeling errors and disturbances are taken into account in the plant. A controller with a simple structure is produced by applying an SVM and L2-gain design. The SVM is used to identify the complicated nonlinear functions and the modeling errors in the plant. The Lagrangian factors in the SVM are obtained using on-line tuning algorithms. L2-gain design is applied to suppress the disturbances. To obtain the optimal parameters in the SVM, then particle swarm optimization (PSO) method is incorporated. The stability and robustness of the close-loop system are confirmed by Lyapunov stability analysis. Numerical simulation is performed to demonstrate the validity of the proposed hybrid controller and its superior performance over a conventional PD controller.

1. Introduction

Ship motion control plays an important role in guaranteeing the safety and economy of a ship in navigation. Well-designed control systems allow ships to sail and perform tasks with adequate reliability and economy in severe seas, which are usually due to uncertain environmental forces induced by wind, waves, and current. Since the ship dynamics are nonlinear essentially and often affected by uncertain environmental forces, it is challenging to achieve a robust and accurate control system for ship motion, which refers to course-keeping, roll stabilization, path following, dynamic positioning, vertical motion control, and station-keeping [1]. Conventional control strategies like PID control can no longer satisfy the requirement of navigation, guidance and control of ships. During recent decades, many advanced control schemes have been developed and successfully applied to ship motion control, including sliding mode variable structure control [2,3]; parameter adaptive control [4,5]; H-infinity robust control [6,7]; neural-network control [8,9]; fuzzy control [10,11]; neuro-fuzzy control [12,13]; line-of-sight based model control [14,15], etc.

To obtain a controller for ship motion, knowledge of the ship’s dynamic characteristics is necessary. This knowledge can be obtained through the mathematical model of ship motion. In ship motion control, usually two kinds of mathematical models are available, i.e., control-design models and high-fidelity models [1]. A precise mathematical model of the plant allows one to design a satisfactory controller by applying some conventional control strategies such as linearization feedback control [16], PID control [17], and adaptive control [18]. However, the ship dynamics is inevitably affected by some uncertainties, which might degrade the controller performance if the uncertainties are not properly compensated for or suppressed in the controller design. Such uncertainties include modeling errors and disturbances. For modeling errors, they refer to parameter errors, ignored high-order modes and unmodelled dynamics of ship motion. Examples of unmodelled dynamics are thrust and torque losses, cross-coupling drag, varying wake, air suction and interaction between the thruster and the hull [19]. For disturbances, they result from random noises in mechanical and electrical equipments, or the environmental forces induced by wind, waves and current. Usually, neither modeling errors nor disturbances can be described by precise mathematical expressions. Therefore, a robust and accurate control scheme is required. H-infinity control is one of the main choices [20,21]. However, there are two difficulties with this method. One is that prior knowledge about the bound of an uncertainty is required, which is challenging sometimes. The other is the difficulty in guaranteeing the control accuracy because of the absence of adaptation to the changing characteristics of uncertainties.

To deal with the uncertainties in a plant, an artificial intelligence (AI)-based controller provides another interesting and effective option. Such an AI-based controller is characterized by nonlinear mapping and on-line learning abilities, which make it possible to identify and compensate for the uncertainties in the plant. Commonly employed AI computing approaches include neural networks (NN), Bayesian probability, fuzzy logic, machine learning, evolutionary computation and genetic algorithms. Combined with AI computation, the controller performance can be greatly improved [22]. In the field of ship motion control, during last decades there have been many NN applications e.g., [8,9,23,24]. Comparatively, SVM is a novel technique of Artificial Intelligence. This kind of learning machine aims to achieve structural risk minimization (SRM) [25] that makes it outperform the algorithms based on empirical risk minimization (ERM) like NN. Better generalization ability can be guaranteed because not only prediction error but also model complexity is considered when performing SRM. Globally optimal solution can also be guaranteed in SVM since convex quadratic programming is performed. Additionally, by applying the kernel trick, the curse of dimensionality can be avoided in SVM, which is usually inevitable in NN. SVM has found wide applications in the field of control engineering and in the past years, there have been some SVM applications to ship motion control. For instance, Liu et al. applied fuzzy SVM to the course-keeping [26]. Jiang et al. proposed SVM based general predictive control (GPC) for ship course-keeping [27]. Luo et al. presented a robust SVM controller for the ship course-keeping [28].

This paper presents a robust controller for the ship course-keeping by using SVM identification and L2-gain design. The investigated plant is composed of the ship dynamics and the dynamics of rudder actuator. Uncertainties are taken account of in the cascaded system, which refer to modeling errors and external disturbance. In the controller design, the modeling errors are identified by SVM while the external disturbance is suppressed by using L2-gain design. To guarantee the robust stability of the controller, a Lyapunov recursive function method is adopted. The main novelty of our approach is the determination of the Lagrangian factors in SVM. Rather than calculate the Lagrangian factors through a matrix equation, one can obtain these factors by solving differential equations proposed. To improve the control accuracy, the parameters in SVM are optimized by incorporating particle swarm optimization (PSO). The rest of the paper is organized as follows. In Section 2, the mathematical model of ship steering, the methodologies of SVM and PSO are described; in Section 3, the details of the controller design and the stability analysis are given; in Section 4, a simulation example is presented to illustrate the validity of the controller designed; and the final section is the concluding remarks.

2. Problem Formulation

2.1. Mathematical Model of Ship Steering

A general six-degree-of-freedom (6DOF) model of ship motion can be described as [29]

$$\begin{array}{l}\dot{\mathit{\eta }}=J(\mathit{\eta })\mathit{v}\\ M_{RB}\dot{\mathit{v}}+C_{RB}(\mathit{v})\mathit{v}=\mathit{\tau }\end{array}\}$$

(1)

where $\mathit{\eta }$ is the generalized position vector with $\mathit{\eta }={[x,y,z,\varphi ,\theta ,\psi ]}^{\mathrm{T}}$, defined in the North-East-Down geographical reference frame; $\mathit{v}$ is the generalized velocity vector with $\mathit{v}={[u,v,w,p,q,r]}^{\mathrm{T}}$, defined in the body-fixed reference frame; $J(\mathit{\eta })$ is the kinematic transformation matrix between the body-fixed reference frame and the North-East-Down geographical reference frame; $M_{RB}$ is the rigid-body mass matrix; $C_{RB}(\mathit{v})$ is the Coriolis-centripetal matrix; $\mathit{\tau }$ denotes a generalized force vector.

System (1) gives a general description of the motion of a marine vessel, a surface ship or an underwater vehicle. In practical applications to ship motion control, this complicated model is often simplified to a model with less than 6DOF. Such a simplified model can be called as a control-design model in which the essential behavior of a dynamic system is captured [1]. More often than not, to study the yaw dynamics and control, a 1DOF model can be adopted. This model is usually named as the response model (sometimes named as Nomoto model or K-T model). Despite simplicity, this model reflects the important response characteristics of the yaw dynamics, from the control input (i.e., the rudder angle) to the system output (i.e., the yaw rate or the heading angle). Practically, this kind of model is commonly preferred in the studies on ship motion control. A first-order nonlinear response model has the form as

$$\begin{array}{l}\dot{\psi }=r\\ T\dot{r}+r+\alpha r^3=K\delta \end{array}\}$$

(2)

where $K,T$ are manoeuvring indices; $\psi$ the heading angle; $r$ the yaw rate; $\alpha$ the coefficient of nonlinear term; $\delta$ the rudder angle. If uncertainties are taken into account, the second subsystem of (2) can be written as

$$T\ddot{\psi }+\dot{\psi }+\alpha {\dot{\psi }}^3+\Delta (\psi ,\dot{\psi })+w=K\delta$$

(3)

where $\Delta (\psi ,\dot{\psi })$ is the modeling error, $w$ is the external disturbance.

In many researches on ship course-keeping control, the actuator dynamics were ignored. It should be noted that if the time constant of the actuator is not much smaller than that of the ship, ignoring the actuator dynamics would result in degradation of the controller performance. Generally, the rudder actuator dynamics can be expressed as

$$T_E\dot{\delta }={\delta }_c-\delta$$

(4)

where $T_E$ is the time constant of rudder; ${\delta }_c$ the demand or execution rudder angle. Combining (3) and (4) through an integrator, a cascaded system is derived as follows, in which ${\delta }_c$ is the system control input; while the rudder angle $\delta$ is an interim state variable and the heading angle $\psi$ is the system output.

$$\begin{array}{l}T\ddot{\psi }+\dot{\psi }+\alpha {\dot{\psi }}^3+\Delta (\psi ,\dot{\psi })+w=K\delta \\ T_E\dot{\delta }={\delta }_c-\delta \end{array}\}$$

(5)

2.2. Support Vector Regression

SVM is a kind of machine learning algorithm for classification and regression. Since it was invented, it has been widely applied in scientific and engineering areas due to its powerful learning abilities. For a MISO system, support vector machines for regression (SVR) have a general form as

$$f(\mathit{x})={\mathit{w}}^{{T}}\cdot \Phi (\mathit{x})+b,\mathit{w}\in {\mathbb{R}}^{\mathit{n}},b\in ℝ$$

(6)

where $\Phi (\cdot )$ is a nonlinear function vector. The above calculation is performed in a so-called high-dimensional feature space to approximate the hidden mapping relationship contained in the original training samples:

$$(\mathit{x}\times y\in {\mathbb{R}}^{\mathit{n}}\times \mathbb{R},) : ({\mathit{x}}_1,y_1),\cdots ,({\mathit{x}}_i,y_i),\cdots ,({\mathit{x}}_l,y_l),i=1,2,\cdots l$$

where $l$ is the number of samples.

To improve feasibility, SVM has been developed during the last decade and some effective modified SVM versions have been proposed, one of which is the Least Squares SVM (LS-SVM). Rather than solve the problem of convex quadratic programming that is required in standard SVM, LS-SVM provides a simpler solution by solving a linear matrix equation [30,31]

$$\left[\begin{array}{cc}0& {\mathit{1}}^{{T}}\\ \mathit{1}& \mathit{Q}+C^{-1}\mathit{I}\end{array}\right]\left[\begin{array}{c}b\\ \mathit{\alpha }\end{array}\right]=\left[\begin{array}{c}0\\ \mathit{y}\end{array}\right]$$

(7)

where

• $\mathit{1}={\left[\begin{array}{ccc}1,& \cdots ,& 1\end{array}\right]}_{1\times l}^{{T}}$;
• $\mathit{I}$ is a $l$-dimensional unit matrix;
• $\mathit{Q}={(q_{ij})}_{l\times l}$ is the kernel matrix with the element $q_{ij}={\Phi }^{{T}}({\mathit{x}}_i)\Phi ({\mathit{x}}_j)$;
• $C$ is the so-called regularization factor, an important hyperparameter that controls the tradeoff between empirical risk and generalization ability of SVM;
• $\mathit{\alpha }={({\alpha }_i)}_{l\times 1}$; $\mathit{y}={(y_i)}_{l\times 1}$.

After the constant $b$ and Lagrangian factors ${\alpha }_i$ are calculated through (7), the estimated function in (6) can be determined by

$$f(\mathit{x})={\displaystyle \sum _{i=1}^{N_s}{\alpha }_{i}K({\mathit{x}}_i,\mathit{x})+b}$$

(8)

where $N_s$ is the number of the support vectors associated with nonzero ${\alpha }_i$. $K({\mathit{x}}_i,{\mathit{x}}_j)$ is the kernel function defined by $K({\mathit{x}}_i,{\mathit{x}}_j)=(\Phi ({\mathit{x}}_i)\cdot \Phi ({\mathit{x}}_j))$.

2.3. Particle Swarm Optimization

Besides SVM, another AI technique, particle swarm optimization (PSO) is also applied in this paper. As a popular parameter optimization tool, PSO is introduced to obtain the optimal parameters in SVM, which exert a vital influence on the performance of SVM. Model selection is an important issue in machine learning, no matter for classifiers or regressors. For SVM, the main concern in model selection is the selection of SVM parameters. Those parameters include insensitivity factor, regularization factor, and kernel parameters. However, how to calculate their optimal values theoretically still remains unsolved. One reason is that it is difficult to determine the VC dimension [32]. A simple and commonly used way to obtain these parameters is the trial-and-error method despite some proposals, for instance analytic parameter selection directly from the training data [33], in-sample and out-of-sample method [34], dynamic particle filter [35], etc. To improve the performance of SVM, this paper makes use of PSO to determine SVM parameters.

The main idea of PSO is to optimize a problem from a population of candidate solutions (so-called particles). These particles are moved towards the optimal position (best solution) through iteration of velocity and position [36]. This optimization method is easy to use since few or no assumptions are required; for example the optimization problem is not required to be differentiable. Moreover, a globally optimal solution can be guaranteed even if the capacity of candidate solutions is very large. In practice, PSO has been proven an effective universal parameter optimizer. There have been some PSO applications to the parameter optimization of SVM, e.g., [37,38,39].

A standard PSO algorithm adopts the following iteration formulae about the particle’s velocity and position

$$V_{id}(k+1)=W\cdot V_{id}(k)+c_1\kappa (p_{id}(k)-X_{id}(k))+c_2\eta (p_{gd}(k)-X_{id}(k))$$

(9)

$$X_{id}(k+1)=X_{id}(k)+V_{id}(k+1),i=1,2,\dots n;d=1,2,\dots q$$

(10)

where $V_{id}(k)$ is the velocity of a particle at the k-th step; n is the number of particles; q is the dimension of particle; $X_{id}$ is the position of a particle; $W$ is the inertia weight; $c_1$ and $c_2$ are acceleration coefficients; $p_{id}$ is the best position of a individual particle while $p_{gd}$ is the swarm or social best position; $\kappa ,\eta$ are both random.

3. Controller Design and Stability Analysis

Usually, the problem of ship course-keeping control can be viewed as a tracking control issue. Without loss of generality, two general assumptions with respect to the continuous cascaded system described by (5) can be given as [A1]. The desired system output signal ${\psi }_d$ is differentiable and the cascaded system is controllable. [A2]. The external disturbance $w(t)$ is bounded with $\left|w(t)\right|\le \overline{w},\forall t\ge 0,\exists \overline{w}>0.$

3.1. Controller Design

First of all, feedback control is applied to the certain terms in the System (3). A desired rudder angle can be defined as

$${\delta }_d=\frac{1}{K}(T{\ddot{\psi }}_d+f(\psi ,\dot{\psi }){\dot{\psi }}_d+u_1)$$

(11)

where $f(\psi ,\dot{\psi })=1+\alpha {\dot{\psi }}^2$, $u_1$ is introduced as an auxiliary controller to deal with the uncertainties in (3).

Three error signals are defined as

$$e={\psi }_d-\psi$$

(12)

$$\xi ={\delta }_d-\delta$$

(13)

$$\zeta =\dot{e}+\lambda e(\lambda >0)$$

(14)

Incorporating (11)–(14) into (3), it follows

$$T\dot{\zeta }=(T\lambda -f(\psi ,\dot{\psi }))\dot{e}-u_1+\Delta (\psi ,\dot{\psi })+w+K\xi$$

(15)

For the actuator dynamics described by (4), the controller can be designed as

$${\delta }_c={\delta }_d+u_2$$

(16)

where $u_2$ is another auxiliary controller. By analogy to (15), it holds

$$T_E\dot{\xi }=T_E{\dot{\delta }}_d-\xi -u_2$$

(17)

Thus, based on (15) and (17), an error system can be obtained as

$$\begin{array}{l}\dot{\zeta }=(\lambda -\frac{f(\psi ,\dot{\psi })}{T})\dot{e}+\frac{\Delta (\psi ,\dot{\psi })}{T}+\frac{w}{T}+\frac{K}{T}\xi -\frac{1}{T}{u}_1\\ \dot{\xi }={\dot{\delta }}_d-\frac{\xi }{T_E}-\frac{1}{T_E}{u}_2\end{array}\}$$

(18)

As can be seen, the above system is linear with respect to the error vector ${\left(\begin{array}{ccc}e& \zeta & \xi \end{array}\right)}^{\mathrm{T}}$. Moreover, the convergence of the error system will depend on the design of two auxiliary controllers, i.e., u₁ and u₂. Lyapunov function method is adopted to obtain proper u₁ and u₂. This method has been widely used in control engineering since it provides an effective and systematic procedure to obtain a stable controller. For the control of a nonlinear complicated system, this method is usually preferred.

For the System (18), a definitely positive Lyapunov function candidate w.r.t. the error vector ${\left(\begin{array}{ccc}e& \zeta & \xi \end{array}\right)}^{\mathrm{T}}$ can be defined as

$$V_1=\frac{1}{2}{e}^2+\frac{1}{2}T{\zeta }^2+\frac{1}{2}{T}_E{\xi }^2$$

(19)

Its derivative w.r.t. time is

$${\dot{V}}_1=-\lambda e^2+\zeta (T\lambda \dot{e}-f(\psi ,\dot{\psi })\dot{e}-u_1+\Delta (\psi ,\dot{\psi })+e+w+K\xi )+\xi (T_E{\dot{\delta }}_d-\xi -u_2)$$

(20)

As can be seen, to obtain the auxiliary controller u₁, one has to deal with the uncertain terms, i.e., $\Delta (\psi ,\dot{\psi })$ and w. For the disturbance $w$, L2-gain design is applied. This approach derives from H-infinity control strategy and is considered as an effective way to suppress uncertain disturbance. Using this method, an L2-gain performance index should be defined as

$${\displaystyle {\int }_0^t{\Vert z\Vert }^2}d\tau \le {\mu }_1^2{\displaystyle {\int }_0^t{w}^{2}d\tau }+{\mu }_2$$

(21)

where ${\mu }_1,{\mu }_2$ are small positive constants, z is introduced as an evaluation signal.

An identity w.r.t. the evaluation signal z and the disturbance w can be obtained as

$${\Vert z\Vert }^2-{\mu }_1^2{w}^2\equiv \frac{1}{4{\mu }_1^2}{\zeta }^2-(\frac{1}{2{\mu }_1}\zeta -{\mu }_{1}w{)}^2+{\Vert z\Vert }^2-\zeta w$$

(22)

Let $z={[r_{1}e,r_2\dot{e}]}^{\mathrm{T}}(r_1,r_2>0)$, the above identity becomes the following inequality

$${\Vert z\Vert }^2-{\mu }_1^2{w}^2\le \frac{1}{4{\mu }_1^2}{\zeta }^2+r_1^2{e}^2+r_2^2({\zeta }^2+{\lambda }^2{e}^2-2\lambda \zeta e)-\zeta w$$

(23)

Combining the above inequality with the derivative (20), one has

$$\begin{array}{r}{\dot{V}}_1+{\Vert z\Vert }^2-{\mu }_1^2{w}^2\le (r_1^2+r_2^2{\lambda }^2-\lambda )e^2+\zeta (T\lambda \dot{e}-f(\psi ,\dot{\psi })\dot{e}+\Delta (\psi ,\dot{\psi })+e\\ +\frac{1}{4{\mu }_1^2}\zeta +r_2^2\zeta -2r_2^2\lambda e-u_1)+\xi (T_E{\dot{\delta }}_d+K\zeta -\xi -u_2)\end{array}$$

(24)

To obtain u₁, one has to account for another uncertainty, i.e., $\Delta (\psi ,\dot{\psi })$. Different from the way used to deal with the disturbance w, the modeling error $\Delta (\psi ,\dot{\psi })$ will be identified and compensated instead of being suppressed by using L2-gain design. This is because although both modeling error and disturbance are uncertainties, the modeling error is usually a continuous signal while the disturbance w might be instantaneous. In this paper SVM is proposed to approximate such a continuous uncertainty. Besides, the term $T_E{\dot{\delta }}_d$ in (24) is also dealt with by SVM. As seen from (24), to obtain u₂, it is necessary to calculate the derivative of δ_d. According to the Equation (11), one has

$${\dot{\delta }}_d=\frac{1}{K}(T{\stackrel{⃛}{\psi }}_d+f(\psi ,\dot{\psi }){\ddot{\psi }}_d+2\alpha \dot{\psi }\ddot{\psi }{\dot{\psi }}_d+{\dot{u}}_1)$$

(25)

The main difficulty is to calculate the derivative of u₁. Therefore, SVM is used to identify this complicated nonlinear function. One can define the two nonlinear functions in (24) and their SVM approximations as

$$f_1=T\lambda \dot{e}-f(\psi ,\dot{\psi })\dot{e}+\Delta (\psi ,\dot{\psi })+e+\frac{1}{4{\mu }_1^2}\zeta +r_2^2\zeta -2r_2^2\lambda e={\displaystyle \sum _{i=1}^{l_1}{\alpha }_{i}K({\mathit{x}}_i,\mathit{x})}+{\eta }_1$$

(26)

$$f_2=T_E{\dot{\delta }}_d+K\zeta -\xi ={\displaystyle \sum _{j=1}^{l_2}{\beta }_{j}K({\mathit{x}}_j,\mathit{x})+{\eta }_2}$$

(27)

where ${\eta }_i(i=1,2)$ is the approximation error and contains the constant b that is defined in (6). Moreover, it is assumed that $\left|{\eta }_1\right|\le {\eta }_{1N},\left|{\eta }_2\right|\le {\eta }_{2N}({\eta }_{1N}>0,{\eta }_{2N}>0)$. ${\alpha }_i$ and ${\beta }_j$ are Lagrangian factors. If the SVM approximates the nonlinear functions f₁ and f₂ very well, the prediction errors ${\eta }_{1,2}$ should be small enough. In that case, the Lagrangian factors, ${\alpha }_i$ and ${\beta }_j$, can be called “ideal” Lagrangian factors. However, these “ideal” factors are actually unable to be obtained in controller design because $\Delta (\psi ,\dot{\psi })$ and ${\dot{\delta }}_d$ in (26) and (27) are unknown uncertainties. Suppose ${\alpha }_i^{\prime }$ and ${\beta }_j^{\prime }$ are the updated approximations of their “ideal” values with high accuracy, the two nonlinear functions f₁ and f₂ can be approximated by

$$f_1\approx {\displaystyle \sum _{i=1}^{l_1}{\alpha }_i^{\prime }K({\mathit{x}}_i,\mathit{x})}$$

(28)

$$f_2\approx {\displaystyle \sum _{j=1}^{l_2}{\beta }_j^{\prime }K({\mathit{x}}_j,\mathit{x})}$$

(29)

Thus, the two auxiliary controllers u₁ and u₂ in (24) can be designed as

$$\begin{array}{c}{u}_1={\displaystyle \sum _{i=1}^{l_1}{\alpha }_i^{\prime }K({\mathit{x}}_i,\mathit{x})}+{\lambda }_{1}T\zeta \\ u_2={\displaystyle \sum _{j=1}^{l_2}{\beta }_j^{\prime }K({\mathit{x}}_j,\mathit{x})}+{\lambda }_2{T}_E\xi \end{array}\}$$

(30)

where ${\lambda }_1$ and ${\lambda }_2$ are positive control gains.

Combining the auxiliary controllers with (11) and (16), the command rudder angle to the cascaded System (5) can be expressed as

$${\delta }_c=\frac{1}{K}(T{\ddot{\psi }}_d+f{\dot{\psi }}_d+{\displaystyle \sum _{i=1}^{l_1}{\alpha }_i^{\prime }K({\mathit{x}}_i,\mathit{x})}+{\lambda }_{1}T\zeta )+{\displaystyle \sum _{j=1}^{l_2}{\beta }_j^{\prime }K({\mathit{x}}_j,\mathit{x})}+{\lambda }_2{T}_E\xi$$

(31)

Figure 1 shows the control framework.

3.2. Stability Analysis

Based on (26), (27) and (30), the inequality (24) becomes

$${\dot{V}}_1+{\Vert z\Vert }^2-{\mu }_1^2{w}^2\le (r_1^2+r_2^2{\lambda }^2-\lambda )e^2+\zeta {\eta }_1+\zeta {\tilde{\mathit{H}}}_1^{\mathrm{T}}{\mathit{\phi }}_1(\mathit{x})-{\lambda }_{1}T{\zeta }^2+\xi {\eta }_2+\xi {\tilde{\mathit{H}}}_2^{\mathrm{T}}{\mathit{\phi }}_2(\mathit{x})-{\lambda }_2{T}_E{\xi }^2$$

(32)

in which

• ${\tilde{\mathit{H}}}_i={\mathit{H}}_i-{\mathit{H}}_i^{\prime },(i=1,2),$
• ${\mathit{H}}_1={\left[\begin{array}{cccc}{\alpha }_1,& {\alpha }_2,& \dots & {\alpha }_{l_1}\end{array}\right]}^{\mathrm{T}}$; ${\mathit{H}}_1^{\prime }={\left[\begin{array}{cccc}{\alpha }_1^{\prime },& {\alpha }_2^{\prime },& \dots & {\alpha }_{l_1}^{\prime }\end{array}\right]}^{\mathrm{T}}$;
• ${\mathit{H}}_2={\left[\begin{array}{cccc}{\beta }_1,& {\beta }_2,& \dots & {\beta }_{l_2}\end{array}\right]}^{\mathrm{T}}$; ${\mathit{H}}_2^{\prime }={\left[\begin{array}{cccc}{\beta }_1^{\prime },& {\beta }_2^{\prime },& \dots & {\beta }_{l_2}^{\prime }\end{array}\right]}^{\mathrm{T}}$;
• ${\mathit{\phi }}_1(\mathit{x})={\left[\begin{array}{cccc}K({\mathit{x}}_1,\mathit{x}),& K({\mathit{x}}_2,\mathit{x}),& \dots & K({\mathit{x}}_{l_1},\mathit{x}),\end{array}\right]}^{\mathrm{T}};$ ${\mathit{\phi }}_2(\mathit{x})={\left[\begin{array}{cccc}K({\mathit{x}}_1,\mathit{x}),& K({\mathit{x}}_2,\mathit{x}),& \dots & K({\mathit{x}}_{l_2},\mathit{x}),\end{array}\right]}^{\mathrm{T}}.$

To make sure the right hand side of (32) be negative, the estimates ${\alpha }_i^{\prime }$ and ${\beta }_i^{\prime }$ are designed by adopting the following on-line tuning algorithms

$$\begin{array}{c}{\dot{\mathit{H}}}_1^{\prime }(t)={\epsilon }_1{\mathit{\phi }}_1(\mathit{x})\zeta (t)-{\epsilon }_2\left|\zeta \right|{\mathit{H}}_1^{\prime }(t)\\ {\dot{\mathit{H}}}_2^{\prime }(t)={\epsilon }_1{\mathit{\phi }}_2(\mathit{x})\xi (t)-{\epsilon }_2\Vert \mathit{\varsigma }\Vert {\mathit{H}}_2^{\prime }(t)\end{array}\}$$

(33)

where ${\epsilon }_1$ and ${\epsilon }_2$ are positive constants. In addition, an augmented error signal is defined as

$$\mathit{\varsigma }={\left[\begin{array}{cc}\zeta & \xi \end{array}\right]}^{{T}}$$

(34)

Based on (19), another Lyapunov function candidate that involves the error vector ${\tilde{\mathit{H}}}_i$ can be defined as

$$V_2=V_1+\frac{1}{2{\epsilon }_1}(tr\{{\tilde{\mathit{H}}}_1^{\mathrm{T}}{\tilde{\mathit{H}}}_1\}+tr\{{\tilde{\mathit{H}}}_2^{\mathrm{T}}{\tilde{\mathit{H}}}_2\})$$

(35)

where $tr\{\cdot \}$ represents the trace of a matrix. Combined with the definition given by (33), the following derivative can be obtained

$$\begin{array}{l}\frac{d}{dt}\left(\frac{1}{2{\epsilon }_1}(tr\{{\tilde{\mathit{H}}}_1^{\mathrm{T}}{\tilde{\mathit{H}}}_1\}+tr\{{\tilde{\mathit{H}}}_2^{\mathrm{T}}{\tilde{\mathit{H}}}_2\})\right)\\ =-\zeta {\tilde{\mathit{H}}}_1^{\mathrm{T}}{\mathit{\phi }}_1-\xi {\tilde{\mathit{H}}}_2^{\mathrm{T}}{\mathit{\phi }}_2+\frac{{\epsilon }_2}{{\epsilon }_1}\left|\zeta \right|tr\{{\tilde{\mathit{H}}}_1^{\mathrm{T}}({\mathit{H}}_1-{\tilde{\mathit{H}}}_1)\}+\frac{{\epsilon }_2}{{\epsilon }_1}\Vert \mathit{\varsigma }\Vert tr\{{\tilde{\mathit{H}}}_2^{\mathrm{T}}({\mathit{H}}_2-{\tilde{\mathit{H}}}_2)\}\end{array}$$

(36)

Let $0\le {\lambda }_3\le 1$ and

$${\lambda }_0=\min \{\begin{array}{ccc}-r_1^2-r_2^2{\lambda }^2+\lambda ,& (1-{\lambda }_3){\lambda }_1,& (1-{\lambda }_3){\lambda }_2,\end{array}\begin{array}{cc}(1-{\lambda }_3){\epsilon }_2\left|\zeta \right|,& (1-{\lambda }_3){\epsilon }_2\Vert \mathit{\varsigma }\Vert \end{array}\}$$

the following inequality w.r.t. (35) holds

$$\begin{array}{cc}\hfill {\dot{V}}_2+{\Vert z\Vert }^2-{\mu }_1^2{w}^2& \le -2{\lambda }_0{V}_2-{\lambda }_3{\lambda }_1{\zeta }^2+\zeta {\eta }_1-{\lambda }_3{\lambda }_2{\xi }^2+\xi {\eta }_2\hfill \\ & +\frac{{\epsilon }_2}{{\epsilon }_1}\left|\zeta \right|tr\{{\tilde{\mathit{H}}}_1^{\mathrm{T}}({\mathit{H}}_1-{\lambda }_3{\tilde{\mathit{H}}}_1)\}+\frac{{\epsilon }_2}{{\epsilon }_1}\Vert \mathit{\varsigma }\Vert tr\{{\tilde{\mathit{H}}}_2^{\mathrm{T}}({\mathit{H}}_2-{\lambda }_3{\tilde{\mathit{H}}}_2)\}\hfill \end{array}$$

(37)

Without loss of generality, the Lagrangian factor vectors H₁ and H₂ can be assumed bounded by ${\Vert {\mathit{H}}_1\Vert }_F\le H_{1M},{\Vert {\mathit{H}}_2\Vert }_F\le H_{2M}(H_{1M},H_{2M}>0)$ where ${\Vert \cdot \Vert }_F$ represents the Frobenius norm with the definition ${\Vert \mathit{A}\Vert }_F=\sqrt{{\displaystyle \sum _{i=1}^m{\displaystyle \sum _{j=1}^n{a}_{ij}^2}}},\forall \mathit{A}={(a_{ij})}_{m\times n}$. Thus, the inequality (37) evolves to

$$\begin{array}{cc}\hfill {\dot{V}}_2+{\Vert z\Vert }^2-{\mu }_1^2{w}^2& \le -2{\lambda }_0{V}_2+\Vert \mathit{\varsigma }\Vert ({\tau }_1+\frac{{\epsilon }_2}{4{\lambda }_3{\epsilon }_1}(H_{1M}^2+H_{2M}^2)\hfill \\ & -{\tau }_2\Vert \mathit{\varsigma }\Vert -\frac{{\lambda }_3{\epsilon }_2}{{\epsilon }_1}{({\Vert {\tilde{\mathit{H}}}_1\Vert }_F-\frac{1}{2{\lambda }_3}{H}_{1M})}^2-\frac{{\lambda }_3{\epsilon }_2}{{\epsilon }_1}{({\Vert {\tilde{\mathit{H}}}_2\Vert }_F-\frac{1}{2{\lambda }_3}{H}_{2M})}^2),(\exists {\tau }_{1,2}>0)\hfill \end{array}$$

(38)

Given a compact set $U_{\mathit{\varsigma }}=\left\{\mathit{\varsigma }|\Vert \mathit{\varsigma }\Vert \le b_{\mathit{\varsigma }},b_{\mathit{\varsigma }}\in {ℝ}^{+}\right\}$, if the tracking error $\mathit{\varsigma }(t)$ falls within the compact set at any time, the above inequality can be reduced to

$${\dot{V}}_2+{\Vert z\Vert }^2-{\mu }_1^2{w}^2\le -2{\lambda }_0{V}_2+b_{\mathit{\varsigma }}({\tau }_1+\frac{{\epsilon }_2}{4{\lambda }_3{\epsilon }_1}(H_{1M}^2+H_{2M}^2))$$

(39)

Let $\epsilon =b_{\mathit{\varsigma }}({\tau }_1+\frac{{\epsilon }_2}{4{\lambda }_3{\epsilon }_1}(W_{1M}^2+W_{2M}^2))$, in the case when the disturbance $w=0$, the above inequality can be simplified to

$${\dot{V}}_2\le -2{\lambda }_0{V}_2+\epsilon$$

(40)

The uniformly ultimately bounded (UUB) stability can be proven [40].

Back to (39), in the case when disturbance exists, i.e., $w\ne 0$, one has

$${\dot{V}}_2+{\Vert z\Vert }^2-{\mu }_1^2{w}^2\le \epsilon$$

(41)

Integrating the above inequality yields

$${\displaystyle {\int }_0^t{\Vert z\Vert }^2}d\tau \le {\mu }_1^2{\displaystyle {\int }_0^t{w}^{2}d\tau }+{\epsilon }_0$$

(42)

where ${\epsilon }_0$ is a small positive constant. Compared with the L2-gain performance index given by (21) it can be said that the disturbance w is suppressed. In other words, the robustness is guaranteed.

If the tracking error $\mathit{\varsigma }(t)$ falls outside the compact set, i.e., $\Vert \mathit{\varsigma }\Vert >b_{\mathit{\varsigma }}$, the following inequality can be obtained

$${\dot{V}}_2+{\Vert z\Vert }^2-{\mu }_1^2{w}^2\le -2{\lambda }_0{V}_2$$

(43)

by appropriately selecting parameters so that

$$b_{\mathit{\varsigma }}\ge ({\tau }_1+\frac{{\epsilon }_2}{4{\lambda }_3{\epsilon }_1}(H_{1M}^2+H_{2M}^2))/{\tau }_2$$

(44)

holds. By analogy to the analysis of (39), in cases either $w=0$ or $w\ne 0$, the stability and robustness of (43) can be guaranteed.

From the above analysis, it can be concluded that the controller designed by (31) and the SVM algorithms given by (33) can guarantee the stability and robustness of the tracking system. It is noted that although only UUB stability can be guaranteed, the control accuracy can be improved by appropriate selection of the parameters in controller, including the control gains and the parameters in SVM. In the study, PSO is incorporated into SVM to obtain optimal SVM parameters.

4. Example Study

To confirm the validity of the controller designed, a numerical simulation is performed combined with a surface ship. Firstly, the state space equation form of (5) can be obtained as

$$\dot{\mathit{X}}=g(\mathit{X})+h(\mathit{X}){\delta }_c$$

(45)

where $\mathit{X}={[\psi \dot{\psi }\delta ]}^{\mathrm{T}}$, $g(\mathit{X})={[g_1{g}_2{g}_3]}^{\mathrm{T}}$, $h(\mathit{X})={[001/T_E]}^{\mathrm{T}}$

$$g_1=\dot{\psi },g_2=-\dot{\psi }-\alpha {\dot{\psi }}^3-\Delta (\psi ,\dot{\psi })-w+K\delta ,g_3=-\delta /T_E$$

The model parameters are [41]

$$\begin{array}{l}T=10.31(\mathrm{s}),\alpha =81.75(\mathrm{s}/{\mathrm{rad}}^2)\\ K=-0.95(1/\mathrm{s}),T_E=0.72(\mathrm{s})\end{array}$$

From the point of view of ship manoeuvrability, this ship is unstable because K is negative. Therefore, an appropriate robust controller is required for the purpose of course-keeping. The controller gains in (31) are selected as ${\lambda }_1={\lambda }_2=2.$ The constant in (14) is selected as $\lambda =0.1$.

For the SVM identification, a Gaussian kernel is adopted as

$$K({\mathit{x}}_{i(j)},\mathit{x})=\exp \left(\frac{-{\Vert {\mathit{x}}_{i(j)}-\mathit{x}\Vert }^2}{2{\sigma }^2}\right)$$

(46)

where $\sigma$ is the kernel width that has an impressive effect on SVM learning ability. The input variable is defined as

$$\mathit{x}=\mathit{x}(t)={\left[\begin{array}{cc}e(t),& \dot{e}\end{array}(t)\right]}^{{T}}$$

(47)

For the support vectors ${\mathit{x}}_{\mathit{i}}$ and ${\mathit{x}}_{\mathit{j}}$, both of them are selected as two successive former input variables which means that $l_1=l_2=2$ for identifiers in (30). Denote $\mathit{x}(t)=\mathit{x}(k)$, the support vectors can be expressed as

$$\begin{array}{l}{\mathit{x}}_1={\left[\begin{array}{cc}e(k-1),& \dot{e}\end{array}(k-1)\right]}^{\mathrm{T}}\\ {\mathit{x}}_2={\left[\begin{array}{cc}e(k-2),& \dot{e}\end{array}(k-2)\right]}^{\mathrm{T}}\end{array}\}$$

(48)

where $k$ is the sampling time.

To obtain a optimized SVM, PSO is incorporated to obtain the optimal $\sigma$ in (46). Moreover, the parameters ${\epsilon }_1$ and ${\epsilon }_2$ in the tuning algorithms of Lagrangian factors given by (33) are also determined by PSO. The fitness in PSO algorithm adopts the root of mean square error (RSME) as

$$\mathrm{RMSE}=\sqrt{\frac{1}{n}{\displaystyle \sum _{i=1}^n{({\dot{\psi }}_{di}-{\dot{\psi }}_i)}^2}}$$

(49)

The number of particles in PSO is set to 20, as is the number of iteration. The initial particle’s positions and velocities are assumed randomly distributed within a given range. Figure 2 presents the results of RMSE, σ, ε₁ and ε₂. As can be seen from the results, the fitness and the parameters converge to stable values, obtained as $\mathrm{RSME}=0.0117,\sigma =0.4356,{\epsilon }_1={\epsilon }_2=0.1$. The three parameters, i.e., σ, ε₁ and ε₂, are used in the SVM based controller.

In the study, the control objective is to track a desired heading angle. Without loss of generality, the desired heading angle ${\psi }_d$ is assumed as a harmonic signal. To verify the stability and robustness of the controller designed, deviation of the heading angle is assumed at the start, which indicates that $\psi (t=0)\ne {\psi }_d(t=0)$. For uncertainties, the modeling error is assumed $\Delta (\psi ,\dot{\psi })=\zeta$ and the external disturbance is assumed as an instantaneous impulse with the amplitude of 1.5 (°/s) at the time $50\mathrm{s}\le t\le 51\mathrm{s}.$ Figure 3 presents the tracking results, including the histories of heading angle $\psi$, yaw rate $r$, demand rudder angle ${\delta }_C$, actual rudder angle $\delta$, and Lagrangian factors $({\alpha }^{\prime },{\beta }^{\prime })$ as well.

From the simulation results, it can be seen that both stability and robustness of the tracking system are achieved. The tracking accuracy is satisfactory as well. Moreover, the control inputs are reasonable and applicable from the point of view of operation.

As a comparison, a conventional PD controller is designed as follows,

$${\delta }_c=k_{1}e+k_2\dot{e}$$

(50)

where k_1,2 are control gains, selected as k₁ = 5, k₂ = 50. Figure 4 presents the simulation results using the PD controller designed. Compared with the performance of the SVM based controller as shown in Figure 3, it can be seen that obvious oscillation happens to the command rudder angle for a PD controller, which implies a limitation for the PD controller from an operational point of view despite the simple structure of the PD controller.

Further comparison is carried out by combining the SVM identification and inverse dynamic compensation that is commonly employed to obtain a feedback controller for a certain system (without uncertainties). By observing the inequality (24), the auxiliary controllers u₁ and u₂ can be designed as

$$\begin{array}{l}{u}_1=T\lambda \dot{e}-f(\psi ,\dot{\psi })\dot{e}+e+\frac{1}{4{\mu }_1^2}\zeta +r_2^2\zeta -2r_2^2\lambda e+{\displaystyle \sum _{i=1}^2{\alpha }_i^{\prime }K({\mathit{x}}_i,\mathit{x})}+{\lambda }_{1}T\zeta ,({\lambda }_1>0)\\ u_2=K\zeta -\xi +{\displaystyle \sum _{i=1}^2{\beta }_i^{\prime }K({\mathit{x}}_i,\mathit{x})}+{\lambda }_2{T}_E\xi ,({\lambda }_2>0)\end{array}\}$$

(51)

which indicates that only the modeling errors $\Delta (\psi ,\dot{\psi })$ and the derivative of desired rudder angle ${\delta }_d$ in the inequality (20) are approximated by SVM, i.e.,

$$\Delta (\psi ,\dot{\psi })={\displaystyle \sum _{i=1}^{l_1}{\alpha }_{i}K({\mathit{x}}_i,\mathit{x})}+{\eta }_1$$

(52)

$$T_E{\dot{\delta }}_d={\displaystyle \sum _{i=1}^{l_2}{\beta }_{i}K({\mathit{x}}_i,\mathit{x})+{\eta }_2}$$

(53)

The terminal control input can be defined as

$$\begin{array}{cc}\hfill {\delta }_c& =(T{\ddot{\psi }}_d+f{\dot{\psi }}_d+T\lambda \dot{e}-f(\psi ,\dot{\psi })\dot{e}+e+\frac{1}{4{\mu }_1^2}\zeta +r_2^2\zeta -2r_2^2\lambda e+{\displaystyle \sum _{i=1}^2{\alpha }_i^{\prime }K({\mathit{x}}_i,\mathit{x})}+{\lambda }_{1}T\zeta )/K\hfill \\ & +K\zeta -\xi +{\displaystyle \sum _{i=1}^2{\beta }_i^{\prime }K({\mathit{x}}_i,\mathit{x})}+{\lambda }_2{T}_E\xi \hfill \end{array}$$

(54)

Comparing controller (54) with (31), one can find that the controller with inverse dynamic compensation in combination with SVM identification (given by (54)) is more complicated than the controller in which the complicated nonlinear function are indentified by SVM (given by (31)). Figure 5 presents the simulation results using controller (54). The parameters are selected as ${\mu }_1=0.25,r_2=1$. In comparison with the results shown in Figure 3, it can be seen that oscillation of the command rudder angle happens at the initial stage and at the time of instantaneous external disturbance despite rapid convergence. Further remarks can be made as follows: (i) no matter whether inverse dynamic compensation is incorporated, the stability and robustness of the control system can be guaranteed due to the use of Lyapunov theory and L2-gain design; (ii) although the response rate of the controller in the case of inverse dynamic compensation (as described by (54)) is better than the case without compensation (as described by (31)), a moderate helming as shown in Figure 3, which can be achieved by using the controller (31), is preferable from the practical point of view. Moreover, as can be seen, the structure of the controller (31) is simpler than the controller (54), which makes sense in practice.

5. Conclusions

In this paper, SVM is applied to the control of a cascaded system with uncertainties with regard to ship course-keeping. The dynamics of the rudder actuator are considered in the plant. To guarantee the real-time performance of the controller, the Lagrangian factors in SVM are determined by on-line tuning algorithms. There are two advantages of the on-line SVM proposed in this paper over conventional on-line SVM. First, no supervised signal is required; second, the capacity of the training samples can be small, which improves the efficiency of the learning machine. Nonlinear functions and an uncertainty (modeling error) in the plant can be effectively identified by the proposed SVM, which helps obtain a controller with simple structure. Another method, L2-gain design, is applied to suppress the disturbance to the plant. The validity of the hybrid controller based on SVM and L2-gain design is confirmed by numerical simulation.

PSO is used to obtain the parameters in SVM. It should be noted that the optimal parameters are actually obtained in off-line way instead of on-line way, which limits its practical application to real-time control. In the next work, an on-line approach to the selection of SVM parameters will be studied.