*Article* **Optimal Damping Concept Implementation for Marine Vessels' Tracking Control**

**Evgeny I. Veremey**

Computer Applications and Systems Department, Saint-Petersburg State University, Universitetskii Prospekt 35, Petergof, 198504 Saint Petersburg, Russia; e\_veremey@mail.ru

**Abstract:** This work presents the results of studies related to the design of stabilizing feedback connections for marine vessels moving along initially given trajectories. As is known, in mathematical formalization, this question leads to a problem of tracking control synthesis for nonlinear and nonautonomous plants. To provide desirable stability and performance features of the closed-loop systems to be synthesized, it is appropriate to use an optimization approach. Unlike the known synthesis methods, which are usually used within the framework of this approach, it is proposed to implement the optimal damping concept first developed by V.I. Zubov in the early 60s of the last century. Modern interpretation of this concept allows constructing numerically effective procedures of control law synthesis taking into account its applicability in a real-time regime. Central attention is focused on the questions connected with practical adaptation of the optimal damping methods for marine control systems. The operability and effectiveness of the proposed approach are illustrated by a practical example of tracking control design.

**Keywords:** marine vessel; tracking controller; stability; functional; optimal damping

#### **1. Introduction**

The nonlinear tracking control of modern marine vessels is one of the most practically significant and theoretically considerable problems in the area of automatic control analysis and design. In particular, tracking control systems are widely used in different branches such as hydrography, inspection of marine constructions, wreck investigation, underwater cable laying, and so on [1,2]. Central theoretical and practical background of tracking control for various moving plants is presented in [3–5] and other fundamental works.

Various issues associated with the design of tracking controllers for marine surface vessels have already been extensively researched and presented in numerous publications (for example, [1,2,6–16]). To evaluate the state of the art in marine tracking control, let us address some modern works presenting this direction of research.

Currently, it is possible to use various ideas to design nonlinear tracking control laws that are reflected in numerous publications, for example, [7,8,13–16]. However, let us note that the mentioned works are not directly oriented to the application of the optimization technique. This makes it difficult to provide the desired dynamic features of the closedloop connections. Now, it seems to be quite evident that the most effective analytical and numerical tool for feedback connections design is the optimization approach. Several aspects of nonlinear tracking control optimization technique are presented in multitudinous scientific publications, including such popular monographs as [4,5,17–21]. As for the simplest stabilization problem, the autopilots with multipurpose structures of optimized control laws are discussed in detail in [22–26].

The sliding mode control technique for marine control applications is discussed in [1,2,6,8]. This direction seems to be quite constructive, but poorly applicable, since it leads to intensive wear of the actuators.

**Citation:** Veremey, E.I. Optimal Damping Concept Implementation for Marine Vessels' Tracking Control. *J. Mar. Sci. Eng.* **2021**, *9*, 45. https:// doi.org/10.3390/jmse9010045

Received: 26 November 2020 Accepted: 29 December 2020 Published: 4 January 2021

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2021 by the author. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

As for the model predictive control (MPC) approach [9], its most significant disadvantage is the large dimension of the minimization problem that is solved at each step of the control process.

Notably, the complexity of this problem is vast because of the many dynamic requirements, restrictions, and conditions that must be satisfied by the chosen control actions.

It should be noted that many scientific works devoted to the tracking control for marine vessels use linear time invariant models of their motion. However, such models are not quite adequate for the problems of deep maneuvering control in angle and positional dynamic variables. Respectively, one of the most important practical difficulties requiring consideration in the design process is the account of nonlinearity and non-autonomy of the control plant model. In most cases, this problem is a source of dynamic instability and poor performance for various systems that were designed based only on linear approximations.

As for the aforementioned optimization approach, its advantages are determined by the flexibility and convenience of modern optimization methods with respect to the relevant practical demands for control design implementation. Certain analytical and numerical methods are used now to compute the optimal controllers for nonlinear and non-autonomous systems subject to various given performance indices. Nevertheless, there is no saying that the optimization approach is recognized overall as a universal instrument to be put into practice for marine tracking controllers design. This can be explained by the presence of some disadvantages connected with computational troubles. Therefore, there exists a vital necessity to develop persistently analytical and numerical methods of control laws design based on optimization ideology adapting to the specific problems for various marine applications.

At present time, numerous approaches are used for a practical solution of these problems [1–12,17,18]. Usually, they are based on Pontryagin's maximum principle, on Bellman's dynamic programming principle, on finite-dimensional approximation in the range of model predictive control (MPC) technique, etc. Unfortunately, all these approaches are connected with the huge extent of calculations that essentially impedes their implementation both for laboratory design activity and real-time regimes of control.

The existence of numerical difficulties motivates us to use other approaches that allow avoiding the aforementioned shortcomings. This work focuses on a different concept that can be applied to design tracking controllers using the theory of optimal damping (OD). This theory, which was first proposed and developed by V.I. Zubov in his works [19–21], provides effective analytical and numerical methods for control calculations with essentially reduced computational consumptions with respect to classical techniques. We believe that this theory was ahead of its time and was undeservedly underutilized for practical control problem solving. This work is one of the attempts to overcome this omission, taking into account the impressive development of modern computer technologies.

In this article, special attention is paid to the control of marine vessels in terms of the forward speed and heading angle. We are considering the regime of the acceleration in order to achieve the specified forward speed with one-time turn along the heading. To achieve desirable stability and performance features of the reference motion, the correspondent tracking controllers were designed based on the OD technique.

The main contribution of this paper is determined by the following statements. First, we propose to use the OD concept to design tracking controllers for marine vessel speed and heading. This has not been the case before. Second, we discuss a new methodology for selecting the functional to be damped, taking into account the desirable features of the closed-loop system in the range of the optimization technique. Hereby, the choice of this functional as the basis is argued by the guarantee asymptotic stability and the desired quality of control processes. Third, we point to the possibility of applying the OD approach to a wide class of nonlinearities in the mathematical model of the vessel. It is noted that this approach can be implemented in real-time regime of a ship's motion. The practical applicability and effectiveness of the proposed technique is illustrated by a controller design for a transport marine ship.

The novelty of the proposed approach with respect to other works lies in the universality and flexibility of proposed nonlinear non-autonomous control laws based on OD computational procedure, which can be implemented in a real-time regime of functioning for marine control plants.

In general, the present study is an extension of the multipurpose approach proposed in [22–27] and developed in [28] with respect to the marine autopilot control laws with the novel structure, taking into account actuators' time delays.

This article is organized as follows. In Section 2, the optimal damping concept for control law synthesis for nonlinear non-autonomous systems is discussed, taking into account certain specific stability and performance requirements for marine control applications. The known background is presented, and the novel ways are proposed to provide tracing controllers synthesis. Section 3.1 is devoted to the OD synthesis problem statement for the forward speed tracking controller and for the tracking autopilot. Central attention is paid to the presentation of mathematical models of the control plant and dynamic requirements for the quality of the closed-loop connection. Section 3.2 presents an exhaustive novel solution for the mentioned synthesis problem based on the optimal damping concept. In Section 3.3, a practical example of tracking controller synthesis is presented to illustrate the applicability and effectiveness of the proposed approach. Finally, Section 4 concludes the article by discussing the overall results of the investigation and indicates how these results can be further developed.

#### **2. Materials and Methods**

As mentioned above, the essence of this paper involves developing an optimal damping technique of tracking control law synthesis for marine vessels with nonlinear and non-autonomous models. In this section, let us first consider the background and some essential features of the OD approach that define the methodological basis of the study.

First of all, let us introduce a commonly used nonlinear robot-like model of the control plant, which represents marine vessel motion for various regimes of its operation [1,2,6]:

$$\begin{cases} \mathbf{M}\dot{\mathbf{v}} + \mathbf{C}(\mathbf{v})\mathbf{v} + \mathbf{D}(\mathbf{v})\mathbf{v} + \mathbf{g}(\eta) = \mathbf{G}\_u \pi + \mathbf{d}, \\ \dot{\eta} = \mathbf{J}(\eta)\mathbf{v}, \end{cases} \tag{1}$$

where vector <sup>ν</sup> ∈ *<sup>R</sup><sup>n</sup>* presents velocities defined in a plant-fixed frame and vector <sup>η</sup> ∈ *<sup>R</sup><sup>n</sup>* contains position dynamical parameters (displacements and angles) in an Earth-fixed frame. External disturbances and controls are presented by the vectors **<sup>d</sup>** ∈ *<sup>R</sup><sup>n</sup>* and <sup>τ</sup> ∈ *<sup>R</sup>m*, respectively. Let us accept that the inertia matrix is positive definite: **M** = **M***<sup>T</sup>* > 0, the matrix of Coriolis-centripetal terms is skew-symmetrical: **<sup>C</sup>**(ν) = <sup>−</sup>**C***T*(ν), and the damping matrix **D**(ν) > 0 is positive definite but non-symmetrical. Vector **g**(η) represents gravitational and buoyancy forces and moments, **J**(η) is the matrix of rotations, and the matrix **G***u* with the constant components reflects controls allocation.

Let us provide a transformation of the body-fixed frame representation (1) to the Earth-fixed one with respect to the vector η. Following [1], this can be done using the following notations:

$$\begin{array}{lcl}\mathbf{M}\_{\eta}(\eta) := \mathbf{J}^{-T}(\eta)\mathbf{M}\mathbf{J}^{-1}(\eta),\\ \mathbf{C}\_{\eta}(\mathbf{v},\eta) := \mathbf{J}^{-T}(\eta)\left[\mathbf{C}(\mathbf{v}) - \mathbf{M}\mathbf{J}^{-1}(\eta)\dot{\mathbf{J}}(\eta)\right]\mathbf{J}^{-1}(\eta),\\ \mathbf{D}\_{\eta}(\mathbf{v},\eta) := \mathbf{J}^{-T}(\eta)\mathbf{D}(\mathbf{v})\mathbf{J}^{-1}(\eta), \; \mathbf{g}\_{\eta}(\eta) := \mathbf{J}^{-T}(\eta)\mathbf{g}(\eta),\\ \mathbf{\tau}\_{\eta} := \mathbf{J}^{-T}(\eta)\mathbf{G}\_{\boldsymbol{\mu}}\mathbf{\tau}, \; \mathbf{d}\_{\eta} = \mathbf{J}^{-T}(\eta)\mathbf{d}.\end{array} \tag{2}$$

In accordance with (2), initial model (1) of the plant takes the form:

$$\ddot{\boldsymbol{\eta}} = -\mathbf{M}\_{\eta}^{-1}(\boldsymbol{\eta}) \Big( (\mathbf{C}\_{\eta}(\mathbf{v}, \boldsymbol{\eta}) + \mathbf{D}\_{\eta}(\mathbf{v}, \boldsymbol{\eta})) \dot{\boldsymbol{\eta}} + \mathbf{g}\_{\eta}(\boldsymbol{\eta}) \ -\boldsymbol{\pi}\_{\eta} - \mathbf{d}\_{\eta} \Big), \ \mathbf{v} = \mathbf{J}^{-1}(\boldsymbol{\eta}) \dot{\boldsymbol{\eta}}.\tag{3}$$

The essence of the tracking control problem is to provide given desirable motion η = η*d*(*t*) of the vessel, using the following state feedback

$$
\mathfrak{r} = \mathfrak{r}(\mathfrak{v}\_{\prime}, \mathfrak{v}, \mathfrak{n}\_{d}(t)),
\tag{4}
$$

which is a nonlinear non-autonomous tracking controller.

Within mathematical formalization, controller (4) must be implemented to provide the zero equilibrium with respect to the tracking error **e**(*t*) := η(*t*) − η*d*(*t*) for the closed-loop system (3), (4), where **d**(*t*) ≡ 0. Naturally, this equilibrium point must be asymptotically stable to guarantee that **e**(*t*) → 0 as *t* → ∞. Let us especially note that the mentioned closed-loop system is nonautonomous, if we have no constant reference motion η*d*(*t*). This gives reasons for us to require the uniform asymptotic stability in global (UGAS) or local (UAS) form. An additional requirement is that the controller (4) provides the desired dynamical features for the closed-loop system (3), (4) under the action of an admissible control τ ∈ *Tu*.

To set the perform of the controller (4) synthesis, we assume that the vector functions <sup>η</sup>*d*(*t*), <sup>ν</sup>(*t*) :<sup>=</sup> **<sup>J</sup>**−1(η*d*(*t*)) **.** η*d*(*t*), and the corresponding τ*d*(*t*) are given. These functions satisfy Equations (1) or (3), i.e., we have:

$$\begin{array}{l} \mathbf{M}\dot{\mathbf{v}}\_{d}(t) + \mathbf{C}(\mathbf{v}\_{d}(t))\mathbf{v}\_{d}(t) + \mathbf{D}(\mathbf{v}\_{d}(t))\mathbf{v}\_{d}(t) + \mathbf{g}(\eta\_{d}(t)) \equiv \mathbf{G}\_{u}\mathbf{r}\_{d}(t),\\ \dot{\mathbf{\eta}}\_{d}(t) \equiv \mathbf{J}(\eta\_{d}(t))\mathbf{v}\_{d}(t). \end{array} \tag{5}$$

Let introduce the following additional notations:

**.**

$$\tilde{\mathbf{x}} := \begin{pmatrix} \tilde{\mathbf{x}}\_1 \\ \mathbf{x}\_2 \end{pmatrix} = \left( \frac{\mathbf{v}}{\eta} \right), \mathbf{f}(\tilde{\mathbf{x}}) := \left( \frac{-\mathbf{M}^{-1}}{\eta} \mathbf{\tilde{L}}^{-1} \underline{\mathbf{C}}(\underline{\mathbf{v}}) + \frac{\mathbf{D}(\mathbf{v})}{\mathbf{J}(\eta)\mathbf{v}} \underline{\mathbf{v}} - \underline{\mathbf{M}}^{-1} \underline{\mathbf{g}}(\underline{\eta}) \right), \mathbf{B} := \left( \frac{\mathbf{M}^{-1} \mathbf{G}\_{\tilde{\mathbf{u}}}}{\mathbf{0}} \underline{\mathbf{v}} \right)$$

which allows us to present Equations (1) and (5) as

$$
\tilde{\mathbf{x}} = \mathbf{f}(\tilde{\mathbf{x}}) + \mathbf{B}\boldsymbol{\tau}\dot{\mathbf{x}}\_d \equiv \mathbf{f}(\mathbf{x}\_d) + \mathbf{B}\boldsymbol{\tau}\_d \tag{6}
$$

supposing that **d**(*t*) ≡ 0.

Let us also consider deflections

$$\mathbf{x} := \tilde{\mathbf{x}} - \mathbf{x}\_d = \begin{pmatrix} \mathbf{x}\_1 \\ \mathbf{x}\_2 \end{pmatrix} := \begin{pmatrix} \mathbf{e}\_\mathbf{v} \\ \mathbf{e} \end{pmatrix} := \begin{pmatrix} \mathbf{e}\_\mathbf{v} \\ \mathbf{e} \end{pmatrix} := \begin{pmatrix} \mathbf{v} - \mathbf{v}\_d \\ \boldsymbol{\eta} - \mathbf{\eta}\_d \end{pmatrix}, \mathbf{u} := \mathbf{r} - \mathbf{r}\_d \tag{7}$$

of the vessel dynamical parameters from the desirable motion.

Then, we can present equations of the vessel in the deflections from the desirable motion. Using notations (7) on the base of (6), we obtain

$$
\dot{\mathbf{x}} = \mathbf{a}(t, \mathbf{x}) + \mathbf{B}\mathbf{u},\tag{8}
$$

where

$$\mathbf{a}(t,\mathbf{x}) = \mathbf{a}(t,\mathbf{e},\mathbf{e}\_{\mathbf{v}}) := \begin{pmatrix} -\mathbf{M}^{-1}[\mathbf{C}(\mathbf{e}\_{\mathbf{v}} + \mathbf{v}\_{d}(t)) + \mathbf{D}(\mathbf{e}\_{\mathbf{v}} + \mathbf{v}\_{d}(t))](\mathbf{e}\_{\mathbf{v}} + \mathbf{v}\_{d}(t)) - \mathbf{M}^{-1}\mathbf{g}(\mathbf{e} + \mathbf{v}\_{d}(t))} \\ -\mathbf{M}^{-1}[\mathbf{C}(\mathbf{e}\_{d}(t)) + \mathbf{D}(\mathbf{v}\_{d}(t))]\mathbf{v}\_{d}(t) - \mathbf{M}^{-1}\mathbf{g}(\mathbf{n}\_{d}(t))} \\ -\mathbf{M}^{-1}[\mathbf{C}(\mathbf{v}\_{d}(t)) + \mathbf{D}(\mathbf{v}\_{d}(t))]\mathbf{v}\_{d}(t) - \mathbf{M}^{-1}\mathbf{g}(\mathbf{n}\_{d}(t))} \end{pmatrix} \tag{9}$$

It is a matter of simple calculations to check that equation (6) has zero equilibrium position, which must be stabilized by the choice of the controller **u** = **u**(*t*, **x**). If this controller is found, then the tracking feedback (4) can be presented as

$$\boldsymbol{\pi} = \boldsymbol{\pi}(t, \mathbf{v}, \boldsymbol{\eta}) = \boldsymbol{\pi}(\boldsymbol{\eta}, \mathbf{v}, \boldsymbol{\eta}\_d) = \mathbf{u}(t, \mathbf{x}) + \boldsymbol{\pi}\_d(t), \; \mathbf{x} := \begin{pmatrix} \mathbf{v} - \mathbf{v}\_d(t) \\ \overline{\mathbf{v}} - \overline{\mathbf{v}\_d(t)} \end{pmatrix}. \tag{10}$$

As for the desirable dynamic features of the closed-loop connection, the most widely used formalized approach is based on minimizing the following integral functional:

$$J = J(t, \mathbf{x}, \mathbf{u}) = \int\_{t\_0}^{\infty} f(\mathbf{r}, \mathbf{x}, \mathbf{u}) \, d\mathbf{r} \,\omega$$

that determines the quality of control processes for a closed-loop system. Here, subintegral function *f* is positive definite, i.e., *f*(*t*, **x**, **u**) ≥ 0 ∀*t* ≥ *t*0, ∀**x**, ∀**u**, *f* = 0 ⇔ **x** = 0, **u** = 0. However, as is well known, there are certain difficulties in directly implementing this approach, which are determined by the huge extent of calculations that essentially impedes their practical implementation.

We propose to overcome the mentioned obstacles using a novel technique based on the OD concept, connected to the OD problem for the synthesis of the control **u**. To state this problem, firstly, let us introduce the functional to be damped as follows:

$$L = L(t, \mathbf{x}, \mathbf{u}) = V(t, \mathbf{x}) + \int\_{t\_0}^{t} f\_0(\mathbf{r}, \mathbf{x}, \mathbf{u}) \, d\mathbf{r},\tag{11}$$

where *V* = *V*(*t*, **x**) is a Lyapunov function candidate, and *f*<sup>0</sup> is a positively defined function. Let us additionally accept that the function *V* satisfies the following conditions:

$$
\alpha\_1(\|\mathbf{x}\|) \le V(t, \mathbf{x}) \le \alpha\_2(\|\mathbf{x}\|) \forall \mathbf{x} \in B\_r \subset E^n, \forall t \ge t\_0. \tag{12}
$$

where α1, α<sup>2</sup> ∈ *K* are comparison functions [4].

The essence of the OD approach consists of the control generation as a function from the current values of variables *t*, **x** in the form

$$\mathbf{u}\_0 = \mathbf{u}\_0(t, \mathbf{x}) = \operatorname\*{argmin}\_{\mathbf{u} \in \mathcal{U}} \mathcal{W}(t, \mathbf{x}, \mathbf{u}) \tag{13}$$

where *<sup>U</sup>* ⊂ *<sup>E</sup><sup>m</sup>* is the metric compact set, and *<sup>W</sup>* is a rate of the functional *<sup>L</sup>* change along the motions of the plant (8):

$$\mathcal{W}(t, \mathbf{x}, \mathbf{u}) := \left. \frac{dL}{dt} \right|\_{\left(8\right)} = \left. \frac{dV}{dt} \right|\_{\left(8\right)} + f\_0(t, \mathbf{x}, \mathbf{u}) = \frac{\partial V}{\partial t}(t, \mathbf{x}) + \frac{\partial V}{\partial \mathbf{x}}(t, \mathbf{x}) \, \mathbf{a}(t, \mathbf{x}) + \frac{\partial V}{\partial \mathbf{x}}(t, \mathbf{x}) \, \mathbf{B} \mathbf{u} + f\_0(t, \mathbf{x}, \mathbf{u}).$$

In other words, it is necessary to find OD controller (13), using an admissible set *<sup>U</sup>* ⊂ *<sup>E</sup><sup>m</sup>* such that ∀**<sup>u</sup>** ∈ *<sup>U</sup>* : <sup>τ</sup>*d*(*t*) + **<sup>u</sup>**(*t*, **<sup>x</sup>**) ∈ *Tu*, ∀*<sup>t</sup>* ≥ 0.


$$\left. \frac{d\mathcal{W}}{d\mathbf{u}} \right|\_{\mathbf{u}=\mathbf{u}\_{0}(t,\mathbf{x})} = 0 \quad \Rightarrow \quad \left[ \frac{\partial V}{\partial \mathbf{x}}(t,\mathbf{x}) \, \mathbf{B} + \frac{\partial F}{\partial \mathbf{u}}(t,\mathbf{x},\mathbf{u}) \right]\_{\mathbf{u}=\mathbf{u}\_{0}(t,\mathbf{x})} = 0. \tag{14}$$

Using the necessary condition (14), one can solve the following nonlinear system

$$\mathbf{a}(t, \mathbf{x}, \mathbf{u}) = \mathbf{b}(t, \mathbf{x})$$

for any point (*t*, **x**) with respect to the vector **u**, where **a**(*t*, **x**, **u**) = col[*ai*(*t*, **x**)], **b**(*t*, **x**) = col[*bi*(*t*, **x**)],

$$a\_i(t, \mathbf{x}, \mathbf{u}) := \left[\frac{\partial F}{\partial \mathbf{u}}(t, \mathbf{x}, \mathbf{u})\right]\_{i} \\ \text{ } b\_i(t, \mathbf{x}) := -\left[\frac{\partial V}{\partial \mathbf{x}}(t, \mathbf{x}) \, \mathbf{B}\right]\_{i} \\ \text{ } i = \overline{1, n}.$$

Based on [4,5], it can be shown that if the function *V* = *V*(*t*, **x**)is such that *W*(*t*, **x**, **u**0(*t*, **x**)) ≤ −α3(**x**) ∀**x** ∈ *Br*, ∀*t* ≥ *t*0, where α<sup>3</sup> ∈ *K*, then the function *V* is control Lyapunov function (CLF) for the plant (8), and the zero equilibrium for the closed-loop system (8), (13) is UAS.

#### **3. Results**

This section is devoted to a practical implementation of the proposed OD approach to nonlinear tracking controllers design for marine vessels. Particular attention is given to the forward speed tracking control law with initially given reference signal.

#### *3.1. Tracking Control Problem for Marine Vessels*

To consider the problems of tracking control for marine applications, let us accept the following widely used [1,2] nonlinear dynamical model of a marine vessel:

$$\begin{array}{ll}\frac{d\overline{V}\_{\mathbf{x}}}{dt} = T\_{\overline{v}}(\overline{V}\_{\mathbf{x}}, \overline{\mathfrak{L}}, \overline{u}\_{1}) + T\_{h}(\overline{V}\_{\mathbf{x}}, \overline{\mathfrak{L}}),\\\frac{d\overline{V}\_{\mathbf{y}}}{dt} = h\_{2}(\overline{V}\_{\mathbf{x}}, \overline{\mathfrak{L}}), \ \frac{d\overline{w}}{dt} = h\_{3}(\overline{V}\_{\mathbf{x}}, \overline{\mathfrak{L}}),\\\frac{d\overline{w}}{dt} = \overline{w}, \ \frac{d\overline{s}}{dt} = \overline{u}\_{2}.\end{array} \tag{15}$$

Here, the following notations are used: *Vx*, *Vy*, and ω are the projections of the velocity vectors on the axes of a vessel-fixed frame *Oxyz* (Figure 1); <sup>ξ</sup> :<sup>=</sup> *Vy* <sup>ω</sup> <sup>ϕ</sup> <sup>δ</sup> *<sup>T</sup>* is the auxiliary vector of dynamical parameters; ϕ is the heading angle, and δ is the vertical rudder deflection. The functions *Tv* and *Th* represent hydrodynamical forces, which are produced by the ship's engine and the water resistance correspondingly.

**Figure 1.** Earth-fixed and vessel-fixed coordinate frames.

The control signals *u*<sup>1</sup> and *u*<sup>2</sup> must be composed by the automatic control system to be designed. The first of them determines a reference surge speed of the vessel, and the second one presents a reference speed of the rudders' deflections.

Assuming that the number of the screw rotations is proportional to a reference surge speed, the function *Tv*(*Vx*, ξ, *u*1) determining a trust force of the screw can be presented in the form

$$T\_v(\nabla\_x, \mathsf{E}, \pi\_1) \equiv \alpha \pi\_1^2 + \beta (\nabla\_x, \mathsf{E}) \pi\_1 + \widetilde{\gamma} (\nabla\_x, \mathsf{E}),\tag{16}$$

We accept here that <sup>α</sup> <sup>≡</sup> *const*, and that the variables <sup>β</sup> and *<sup>γ</sup>*<sup>1</sup> are treated as known functions of the dynamical parameters *Vx*, ξ.

Let us introduce new vector variables **<sup>h</sup>**(*Vx*, <sup>ξ</sup>) :<sup>=</sup> *<sup>h</sup>*2(*Vx*, <sup>ξ</sup>) *<sup>h</sup>*1(*Vx*, <sup>ξ</sup>) <sup>ξ</sup><sup>3</sup> <sup>0</sup> *<sup>T</sup>* , **b** = 0001 *<sup>T</sup>* , and, taking into account (16), we can rewrite nonlinear model (1) of the vessel dynamic as follows:

$$\begin{array}{l}\frac{d\nabla\_{\mathbf{x}}}{dt} = \alpha\overline{u}\_{1}^{2} + \beta(\nabla\_{\mathbf{x}}, \overline{\mathbf{\xi}})\overline{u}\_{1} + \gamma(\overline{\nabla}\_{\mathbf{x}}, \overline{\mathbf{\xi}}),\\\frac{d\overline{\mathbf{\xi}}}{dt} = \mathbf{h}(\nabla\_{\mathbf{x}}, \overline{\mathbf{\xi}}) + \mathbf{b}\overline{u}\_{2}.\end{array}\tag{17}$$

Now we specify a certain reference motion *Vx*ρ(*t*), ξρ(*t*), *u*1ρ(*t*), and *u*2ρ(*t*) of the plant (17), satisfying the identities

$$\begin{array}{l}\frac{dV\_{x\rho}}{dt} \equiv \mathfrak{a}u\_{1\rho}^{2} + \mathfrak{k}(V\_{x\rho\prime}\mathfrak{k}\_{\rho})u\_{1\rho} + \gamma(V\_{x\rho\prime}\mathfrak{k}\_{\rho}),\\\frac{d\underline{k}\_{\rho}}{dt} \equiv \mathfrak{h}(V\_{x\rho\prime}\mathfrak{k}\_{\rho}) + \mathfrak{b}u\_{2\rho}.\end{array} \tag{18}$$

Using systems (17) and (18), we can present equations of vessel dynamic in deflections from the desirable reference motion of the form

$$\begin{split} \frac{dV\_x}{dt} + \frac{dV\_{\rm{sp}}}{dt} &= \alpha \left( u\_1 + u\_{1\rho} \right)^2 + \beta \left( V\_x + V\_{x\rho}, \xi + \xi\_\rho \right) \left( u\_1 + u\_{1\rho} \right) + \\ &+ \gamma \left( V\_x + V\_{x\rho}, \xi + \xi\_\rho \right), \\ \frac{d\xi}{dt} + \frac{d\xi\_\rho}{dt} &= \mathbf{h} \left( V\_x + V\_{x\rho}, \xi + \xi\_\rho \right) + \mathbf{b} \left( u\_2 + u\_{2\rho} \right), \end{split} \tag{19}$$

where

$$V\_{\mathbf{x}} := \overline{V}\_{\mathbf{x}} - V\_{\mathbf{x}\boldsymbol{\rho}\_{\prime}} \ \xi := \overline{\xi} - \xi\_{\boldsymbol{\rho}\_{\prime}} \ u\_{1} := \overline{u}\_{1} - u\_{1\boldsymbol{\rho}\_{\prime}} \ u\_{2} := \overline{u}\_{2} - u\_{2\boldsymbol{\rho}\_{\prime}}.\tag{20}$$

One can see that *Vx*<sup>ρ</sup> = *Vx*ρ(*t*), ξ<sup>ρ</sup> = ξρ(*t*), *u*1<sup>ρ</sup> = *u*1ρ(*t*), and *u*2<sup>ρ</sup> = *u*2ρ(*t*) are known functions of *t*: using new correspondent notations β1, *γ*1, **h**1, we can rewrite (19) as follows:

$$\begin{array}{l} \frac{dV\_{\mathbf{x}}}{dt} = \alpha u\_1^2 + \beta\_1 (t, V\_{\mathbf{x}\prime} \mathfrak{E}) u\_1 + \gamma\_1 (t, V\_{\mathbf{x}\prime} \mathfrak{E}),\\ \frac{d\mathfrak{E}}{dt} = \mathbf{h}\_1 (t, V\_{\mathbf{x}\prime} \mathfrak{E}) + \mathbf{b} u\_2 \end{array} \tag{21}$$

stating that the resulting system (21) has a zero-equilibrium position.

Let us especially note that, unlike (17), this system is not only nonlinear, but also non-autonomous. This is due to the explicit introduction of the time-dependent reference signals into the vessel dynamics equations.

The purpose of the control design procedure is to construct the following stabilizing feedback controllers in deviations:

$$
\mu\_1 = \mu\_1(t, V\_\mathbf{x}, \mathfrak{E}), \; \mu\_2 = \mu\_2(t, V\_\mathbf{x}, \mathfrak{E}), \tag{22}
$$

such that the zero equilibrium of the closed-loop connection (21), (22) is asymptotically stable. This means that if the motion of the initial plant (15) takes place under the action of tracking controllers of the form

$$
\overline{u}\_1 = u\_{1\rho}(t) + u\_1(t, V\_\mathbf{x}, \mathfrak{E}), \\
\overline{u}\_2 = u\_{2\rho}(t) + u\_2(t, V\_\mathbf{x}, \mathfrak{E}), \tag{23}
$$

starting in some neighborhood of the point {*Vx*ρ(0), ξρ(0)}, then this motion tends to the reference one as *t* becomes infinite.

As for the performance of control processes, they are usually formalized mathematically using certain functionals, which are given on the motions of the closed-loop system (21), (22). Currently, the commonly used approach to design stabilizing controllers (22) is setting and solving the following optimization problem

$$J = J(u\_1, u\_2) \to \min\_{\{u\_1, u\_2\} \in \Omega} \tag{24}$$

based on the integral functional

$$J = J(\mu\_1, \mu\_2) = \int\_{t\_0}^{\infty} f^\*(t, V\_x(t), \xi(t), \mu\_1(t), \mu\_2(t)) \, dt,\tag{25}$$

where *U*1 is the set of admissible pairs {*u*1, *u*2} and subintegral function *f* <sup>∗</sup> is positive definite for all its arguments.

In contrast to the problems (24), (25) with traditional methods of its solving, as mentioned above, we propose to use novel approach based on the OD concept.

Let us especially note that both the solution of the problem (24) and the solution of the OD problem significantly depend on the initially given mathematical model (21) of a marine vessel. Naturally, this is evidence that this model cannot be formed accurately, which raises very important questions about the robust features of the closed-loop system to be designed. This problem is quite solvable for various types of marine vessels: the paper [29] is an example. However, the study of the robust features of tracking OD controllers presents an independent problem to be addressed in future studies.

#### *3.2. Forward Speed OD Tracking Controller Design*

In general, the synthesis of the two stabilizing controllers (22) using OD technique can be carried out simultaneously. Nevertheless, one can also apply a combined approach to the feedback design for the considered plant (21) with two control channels. In the range of this approach, the first control is formed based on the OD problem, and the second one can be designed in any other way, providing desirable performance features. However, a joint closed-loop system must have zero equilibrium, and this motion must be asymptotically stable.

Realizing this idea, let us accept the dynamic controller for rudders (marine autopilot) in the following form:

$$\begin{array}{l} \frac{d\mathbf{z}}{dt} = \mathbf{y}\_z(t, \mathbf{z}, V\_{\mathbf{x}}, \mathbf{\xi}),\\ \mu\_2 = \mathbf{g}\_z(t, \mathbf{z}, V\_{\mathbf{x}}, \mathbf{\xi}), \end{array} \tag{26}$$

where **<sup>z</sup>** <sup>∈</sup> *<sup>E</sup>*<sup>ν</sup> is the state vector of this controller.

In particular, the feedback (26) may have a multipurpose structure, which is presented in detail in [22,27] for linear time-invariant case. Its implementation taking into account control time-delay is investigated in the paper [28].

The equations of the control plant now take the form

$$\begin{array}{l}\frac{dV\_{\mathbf{r}}}{dt} = \alpha u\_{1}^{2} + \beta\_{1}(t, V\_{\mathbf{x}}, \xi)u\_{1} + \gamma\_{1}(t, V\_{\mathbf{x}}, \xi),\\\frac{d\xi}{dt} = \mathbf{h}\_{1}(t, V\_{\mathbf{x}}, \xi) + \mathbf{b}\varrho\_{z}(t, \mathbf{z}, V\_{\mathbf{x}}, \xi),\\\frac{d\mathbf{z}}{dt} = \mathbf{y}\_{z}(t, \mathbf{z}, V\_{\mathbf{x}}, \xi).\end{array} \tag{27}$$

In order to synthesize a feedback for the first control channel, i.e., design the OD forward speed stabilizing controller, let us introduce the functional to be damped as

$$L = L(t, V\_{\mathbf{x}}, \xi, \mathbf{z}, u\_1) = V(V\_{\mathbf{x}}, \xi, \mathbf{z}) + \int\_{t\_0}^{t} \lambda\_1^2 u\_1^2 \, d\mathbf{\tau}.\tag{28}$$

Let us take as a Lyapunov function candidate the following sum of quadratic forms:

$$V = V(V\_{\mathbf{x}}, \xi, \mathbf{z}) = \lambda^2 V\_{\mathbf{x}}^2 + \xi^T \mathbf{Q} \xi + \mathbf{z}^T \mathbf{Q}\_1 \mathbf{z}, \ \mathbf{Q} > 0, \ \mathbf{Q}\_1 > 0. \tag{29}$$

Based on (28), (29), we can state the correspondent OD problem of the form

$$\mathcal{W}(t, V\_{\mathbf{x}}, \xi, \mathbf{z}, u\_1) \to \min\_{u\_1 \in \mathcal{U}\_1} \tag{30}$$

where the rate *W* of the functional *L* change is determined as follows:

$$\begin{split} \mathcal{W} &= \mathcal{W}(t, V\_{\mathbf{x}}, \boldsymbol{\xi}, \mathbf{z}, \boldsymbol{u}\_{1}) := \left. \frac{d\mathcal{L}}{dt} \right|\_{\begin{pmatrix} 2\varPi \end{pmatrix}} = \left. \frac{dV}{dt} \right|\_{\begin{pmatrix} 2\varGamma \end{pmatrix}} + \lambda\_{1}^{2} \boldsymbol{u}\_{1}^{2} = \left. \frac{\partial V}{\partial V\_{\mathbf{x}}} \frac{dV\_{\mathbf{x}}}{dt} + \frac{\partial V}{\partial \mathbf{k}} \frac{d\mathbf{k}}{dt} + \frac{\partial V}{\partial \mathbf{z}} \frac{d\mathbf{z}}{dt} + \lambda\_{1}^{2} \boldsymbol{u}\_{1}^{2} \\ &= 2\lambda^{2} V\_{\mathbf{x}} \Big[ \alpha \mathbf{u}\_{1}^{2} + \beta\_{1} (t, V\_{\mathbf{x}}, \boldsymbol{\xi}) \boldsymbol{u}\_{1} + \gamma\_{1} (t, V\_{\mathbf{x}}, \boldsymbol{\xi}) \Big] + 2\xi^{T} \mathbf{Q} \frac{d\mathcal{L}}{dt} + 2\mathbf{z}^{T} \mathbf{Q}\_{1} \frac{d\mathbf{z}}{dt} + \lambda\_{1}^{2} \boldsymbol{u}\_{1}^{2}. \end{split} \tag{31}$$

Let us solve the problem (30), taking into account (31) and assuming that the extreme is achieved at the inner point of the set *U*1: with necessity we obtain

$$\frac{dW}{du\_1} = 2\lambda^2 V\_\mathbf{x} (2\pi u\_1 + \beta\_1 (t, V\_\mathbf{x}, \xi)) + 2\lambda\_1^2 u\_1 = 0$$

that determines the following controller

$$u\_1 = u\_{10}(t, V\_{\mathbf{x}'} \mathfrak{E}) = -\frac{\lambda^2 V\_{\mathbf{x}} \mathfrak{E}\_1(t, V\_{\mathbf{x}'} \mathfrak{E})}{2\lambda^2 V\_{\mathbf{x}} \mathfrak{a} + \lambda\_1^2}.\tag{32}$$

Since, according to formulae (15)–(21), we have

$$\begin{aligned} \beta\_1 &= \beta\_1(t, V\_\ge, \xi) := 2\alpha u\_{1\rho}(t) + \beta \left(V\_\ge + V\_{\ge \rho}(t), \xi + \xi\_\rho(t)\right), \\ \beta(V\_\ge + V\_{\ge \rho}(t), \xi + \xi\_\rho(t)) &= \beta\_0 \sqrt{\left(V\_\ge + V\_{\ge \rho}(t)\right)^2 + \left(V\_\ge + V\_{y\rho}(t)\right)^2}, \ \beta\_0 = \text{const} \end{aligned}$$

we arrive from (32) at the following OD controller:

$$u\_{10}(t, V\_x, \xi) = -\frac{\lambda^2 V\_x \left(2\alpha u\_{1\rho}(t) + \beta\_0 \sqrt{\left(V\_x + V\_{x\rho}(t)\right)^2 + \left(V\_y + V\_{y\rho}(t)\right)^2}\right)}{2\lambda^2 V\_x \alpha + \lambda\_1^2},\tag{33}$$

where *Vy* = 1000 ξ, *Vy*ρ(*t*) = 1000 ξρ(*t*).

It is necessary to note that if the value *u*10(*t*, *Vx*, ξ) is out of an inner part of the admissible set *U*1, we must use the control signal *u*<sup>1</sup> = 0.01*P*ρ*u*1ρ(*t*) · sign(*u*10(*t*, *Vx*, ξ)) instead of (33).

It is necessary to note that practical problems involve situations where the right part of the equation for forward speed has a more complex structure than for the system (17). In general, this equation can be presented as

$$\frac{d\nabla\_{\mathbf{x}}}{dt} = F\_{\mathbf{x}}(\overline{V}\_{\mathbf{x}}, \overline{\mathbf{\xi}}\_{\mathbf{\prime}} \overline{u}\_{1}) .$$

which results in the correspondent equation for the system (21) in deflections:

$$\frac{dV\_{\mathbf{x}}}{dt} = G\_{\mathbf{x}}(t\_{\prime}V\_{\mathbf{x}\prime}\boldsymbol{\xi}\_{\prime}\boldsymbol{u}\_{1})\_{\prime\prime}$$

where *Gx*(*t*, *Vx*, ξ, *u*1) := *Fx Vx* + *Vx*ρ(*t*), ξ + ξρ(*t*), *u*<sup>1</sup> + *u*1ρ(*t*) − *Fx Vx*ρ(*t*), ξρ(*t*), *u*1ρ(*t*) . In this case, instead of (31) we have

$$\mathcal{W}(t, V\_{\mathbf{x}}, \xi, \mathbf{z}, \boldsymbol{\mu}\_{1}) = 2\lambda^{2} V\_{\mathbf{x}} G\_{\mathbf{x}}(t, V\_{\mathbf{x}}, \xi, \boldsymbol{\mu}\_{1}) + 2\xi^{T} \mathbf{Q} \frac{d\xi}{dt} + 2\mathbf{z}^{T} \mathbf{Q}\_{1} \frac{d\mathbf{z}}{dt} + \lambda\_{1}^{2} \boldsymbol{\mu}\_{1}^{2}.$$

and the analytical search for stationary points of the function *W* becomes problematic. Nevertheless, it is always possible to consider a question about the numerical solution of OD problem (30) for each aggregate of fixed parameters *t*, *Vx*, ξ. Moreover, it is always possible to pose a finite dimensional minimization problem

$$u\_1 = u\_{10}^n(t, V\_{\ge}, \xi) = \arg\min\_{u\_1 \in \mathcal{U}\_{1n}} \left( \mathcal{g}\_1(t, V\_{\ge}, \xi, u\_1) + \lambda\_1^2 u\_1^2 \right),\tag{34}$$

on the finite net *U*1*<sup>n</sup>* ⊂ *U*<sup>1</sup> for any compact set *U*1. This problem should be solved at the time *t* for the fixed parameters *t*, *Vx*(*t*), ξ(*t*). Obviously, with a sufficiently large number of points for the set *U*1*n*, the solution of the problem (34) tends to solution for the same problem on the set *U*1.

Finally, in accordance with formula (23), we can compose the tracking controller of the form

$$
\overline{u}\_1 = \overline{u}\_{01}(t, V\_{\ge \prime} \mathfrak{E}) := u\_{1\rho}(t) + u\_{10}(t, V\_{\ge \prime} \mathfrak{E})\_{\prime},
$$

where the stabilizing part *u*10(*t*, *Vx*, ξ) is determined by equality (33). One can see that the resulting representation for OD forward speed tracking controller is as follows:

$$\overline{u}\_{1} = \overline{u}\_{01}(t, \overline{V}\_{x}, \overline{V}\_{y}) = \frac{\lambda\_{1}^{2} u\_{1\rho}(t) - \lambda^{2} (\overline{V}\_{x} - V\_{x\rho}(t)) \left\| \beta\_{0} \sqrt{\overline{V}\_{x}^{2} + \overline{V}\_{y}^{2}} \right. \tag{35}$$

The current values of the dynamic parameters *Vx*(*t*), *Vy*(*t*) must be measured to implement the controller (35).

#### *3.3. Numerical Example of OD Tracking Controller Synthesis*

To illustrate a practical implementation of the proposed OD approach, let us consider a practical example of forward speed tracking controller design for the transport ship with a displacement of 3500 tons, a length of 110 m, a width of 14 m, and an immersion of 5 m. As a mathematical model of the plant, let us accept equations (17), which are presented in [1,2] with the given functions β(*Vx*, ξ), *γ*(*Vx*, ξ), **h**(*Vx*, ξ) and parameter α.

First, we define the reference motion of the ship as a partial solution of system (17) with the initial conditions *Vx*(0) = 2 *m*/*s*, ξ(0) = 0. Let us accept the reference forward speed control as follows:

$$u\_{1\rho}(t) = \begin{cases} 2 + 0.08t, & \text{if } 0 \le t \le 100s, \\ 10, & \text{if } t \ge 100s. \end{cases} \tag{36}$$

To determine the reference control for the rudders, let us initially construct a model feedback as a controller

$$u\_m = k\_1 \mathbf{x}\_{m1} + k\_2 \mathbf{x}\_{m2} + k\_3 (\mathbf{x}\_{m3} - \boldsymbol{\varphi}\_z) + k\_4 \mathbf{x}\_{m4} \\ := \mathbf{k}\_m \mathbf{x}\_m - k\_3 \boldsymbol{\varphi}\_z. \tag{37}$$

fulfilling the command-heading signal ϕ*z*. The basis item **k***m***x***<sup>m</sup>* in (37) stabilizes Linear Time Invariant (LTI) plant

$$\frac{d\mathbf{x}\_m}{dt} = \mathbf{A}\_m \mathbf{x}\_m + \mathbf{b}\_m u\_m \tag{38}$$

which is a result of the plant (17) linearization in the neighborhood of the origin for a fixed forward speed *Vx*. Let us especially notice that the linear model (38) is used only for the synthesis of controller (37) but is not used for simulation and for performance indices computation: controller (37) closes full initial model (17) of the vessel.

Accepting *Vx* = 10 m/s , we obtain

$$\mathbf{A}\_{m} = \begin{pmatrix} a\_{11} & a\_{12} & 0 & b\_{1} \\ a\_{21} & a\_{22} & 0 & b\_{2} \\ 0 & 1 & 0 & 0 \\ 0 & 0 & 0 & 0 \end{pmatrix}, \ \mathbf{b}\_{m} = \begin{pmatrix} 0 \\ 0 \\ 0 \\ 1 \end{pmatrix} \\ a\_{11} = -0.0936, \quad a\_{12} = -6.34, \quad b\_{1} = -0.190, \\ b\_{21} = -0.00480, \quad a\_{22} = -0.717, \quad b\_{2} = 0.0160.$$

We design the basis control *um* = **k***m***x***<sup>m</sup>* for plant (38) with aforementioned parameters as the Linear Quadratic Regulator (LQR) optimal controller with respect to the functional

$$J\_m = J(u\_m) = \int\_0^\infty \left( \mathbf{x}\_m^T \mathbf{Q}\_m \mathbf{x}\_m + \lambda\_m u\_m^2 \right) dt\_\prime \tag{39}$$

where **<sup>Q</sup>***<sup>m</sup>* = diag 0 1.975 0.0250 0 , *λ<sup>m</sup>* = 60. Let us especially remark that the choice of the presented parameters for the functional (39) is determined by the initially given requirements to a dynamic quality of the control processes for the nonlinear time varying system. These requirements provide desirable settling time and overshoot for the closed-loop connection.

As a result of computations, we get the constant vector **k***<sup>m</sup>* = *k*<sup>1</sup> *k*<sup>2</sup> *k*<sup>3</sup> *k*<sup>4</sup> , where

$$k\_1 = 0.00234, \ k\_2 = -0.0495, \ k\_3 = -0.0204, \ k\_4 = -0.0497.$$

Let us substitute the reference heading control *u*1ρ(*t*) and the rudders feedback control

$$
\overline{u}\_2 = k\_1 \overline{\xi}\_1 + k\_2 \overline{\xi}\_2 + k\_3 (\overline{\xi}\_3 - \varphi\_z) + k\_4 \overline{\xi}\_4,\\
\varphi\_z = 90^\circ,\tag{40}
$$

into Equation (17). After integration, we obtain the functions *u*2<sup>ρ</sup> = *u*2ρ(*t*), *Vx*<sup>ρ</sup> = *Vx*ρ(*t*), and ξ<sup>ρ</sup> = ξρ(*t*) = *Vy*ρ(*t*) ωρ(*t*) ϕρ(*t*) δρ(*t*) *T* , which determine given reference motion of the vessel. The graphs of aforementioned functions are presented in Figures 2 and 3.

**Figure 2.** Graphs of the reference functions *u*1ρ(*t*), *u*2ρ(*t*), *Vx*ρ(*t*), and δρ(*t*).

**Figure 3.** Graphs of the reference functions ϕρ(*t*), *Vy*ρ(*t*), and ωρ(*t*).

To provide simulation process, instead of (23), (26), we accept the following simple feedback control law for a marine autopilot:

$$
\overline{u}\_2 = k\_1 V\_y + k\_2 \overline{w} + k\_3 (\overline{\varphi} - \varphi\_z) + k\_4 \overline{\delta}, \tag{41}
$$

where ϕ*<sup>z</sup>* is the command-heading signal. This controller, corresponding to (40), can be directly used for actuators of the vertical rudders.

Thus, one can see that, as a result, initial control plant (17) is closed by the tracking controllers (35) and (41). The current values of the dynamic parameters *Vx*(*t*), *Vy*(*t*), ω(*t*), ϕ(*t*), and δ(*t*) must be measured by the corresponding sensors for the actual implementation of these controllers.

For simulation of the closed-loop system dynamics, the following parameters values are accepted: *λ*<sup>2</sup> = 100, *λ*<sup>2</sup> <sup>1</sup> = 1, α = 0.00462, β<sup>0</sup> = −0.00322. In addition, let us take into account the restrictions δ(*t*) ≤ *dm* = 35◦ and |*u*2(*t*)| ≤ *um* = 3 ◦/c.

To illustrate the practical applicability of the proposed approach, let us simulate the control processes for the closed-loop connection. The aim is to make the proposed approach comparable to other methods. This determined the choice of design parameters and regimes of vessel's motion. These regimes represent the most popular options for movement on quiet water and under the action of sea waves.

The results of simulation are presented in Figures 4–6 as the graphs of corresponding functions, which reflect control signals and the vessel's state variables for the transient process. This process is determined by the aforementioned reference motion, which is realized with the help of the designed tracking controllers. Initial conditions for all variables are zero with the exception of forward speed and heading angle. By these variables, the initial conditions *Vx*(0) = 4 m/s and ϕ(0) = −10◦ are accepted to distinguish them from the reference motion.

**Figure 4.** First control signals *u*<sup>1</sup> = *u*1*b*(*t*), *u*1ρ(*t*), and *u*1(*t*) = *u*1(*t*) − *u*1ρ(*t*).

**Figure 5.** Forward speed *Vx* = *Vxb*(*t*), *Vx*ρ(*t*), and *Vx*(*t*) = *Vx*(*t*) − *Vx*ρ(*t*).

**Figure 6.** Rudders deflections δ = δ*b*(*t*), δρ(*t*), and δ(*t*) = δ(*t*) − δρ(*t*).

Let us note that the dynamical quality of the presented transient process seems to be quite satisfactory. In addition, it is suitable to illustrate the dynamics of the closed-loop connection under the action of sea wave external disturbance. Figure 7 shows the same control process for the forward speed *Vx*, taking into account presence the approximate representation of sea waves with an intensity of 5 on the Beaufort scale.

**Figure 7.** Forward speed *Vx* = *Vxb*(*t*), *Vx*ρ(*t*), and *Vx*(*t*) = *Vx*(*t*) − *Vx*ρ(*t*) under sea waves action.

Finally, let us analyze the stability of the closed-loop system (17), (35), (41) with synthesized tracking controllers. One can easily see that in order to provide aforementioned analysis, it is sufficient to consider the zero-equilibrium stability for the closed-loop system (27), (33) presented in deflections from the reference motion.

It is obvious that systems (27), (33) have zero equilibrium, and to investigate its stability features, let us introduces the following Lyapunov function candidate:

$$V = V(\mathbf{x}) = \lambda^2 V\_x^2 + \xi^T \mathbf{S} \boldsymbol{\xi} = \mathbf{x}^T \mathbf{S}\_1 \mathbf{x}, \; \mathbf{S}\_1 = \left(\cdot, \frac{\lambda^2}{0}, \frac{\lambda}{\cdot} \frac{0}{\mathbf{S}} \cdot \right), \; \mathbf{x} := \begin{pmatrix} V\_x \\ \xi \end{pmatrix}. \tag{42}$$

Here, the symmetrical matrix **S** > 0 is a solution of the algebraic Riccati equation, which is used in the range of LQR controller (37) synthesis:

$$\mathbf{S} = \begin{pmatrix} 0.0107 & -0.243 & -0.102 & -0.141 \\ -0.243 & 7.04 & 2.59 & 2.97 \\ -0.102 & 2.59 & 1.27 & 1.22 \\ -0.141 & 2.97 & 1.22 & 2.98 \end{pmatrix}.$$

Let us notice that the introduced function *V*(**x**) satisfy the relationships

$$
\alpha\_1(\|\mathbf{x}\|) \le V(\mathbf{x}) \le \alpha\_2(\|\mathbf{x}\|), \forall \mathbf{x} \in E^5,\tag{43}
$$

with the following *K*-class functions:

$$\alpha\_1(\left\|\mathbf{x}\right\|) = \sigma\_{\min}(\mathbf{S}\_1) \cdot \left\|\mathbf{x}\right\|^2, \ \alpha\_2(\left\|\mathbf{x}\right\|) = \sigma\_{\max}(\mathbf{S}\_1) \cdot \left\|\mathbf{x}\right\|^2.$$

where *σ*min = 0.00115 and *σ*max = 100 are the minimum and maximum eigenvalues of the matrix **S**1, respectively.

Using the additional function α3(**x**) ≡ **x** ∈ *K*, one can check if the following inequality is valid

$$\mathcal{W}(t, \mathbf{x}) = \frac{dV}{dt}\Big|\_{(2\mathcal{T}),(3\mathbf{3})} \le -\alpha\_{\mathbb{B}}(\|\mathbf{x}\|), \forall t \ge 0, \ \forall \mathbf{x} \in B\_{\mathbf{x}}.\tag{44}$$

where *Bx* is a box, determined by the relationships:

$$B\_{\mathbf{x}}(t) = \left\{ \mathbf{x} \in E^5 : |\mathbf{x}\_i| \le \mathbf{x}\_{\text{im}}(t), \ i = \overline{1,5} \right\},$$

*x*1*m*(*t*) = 0.01*P*ρ|*Vx*ρ(*t*)|, *x*2*m*(*t*) = 0.01*P*<sup>ρ</sup> *Vy*ρ(*t*) ,*x*3*m*(*t*) = 0.01*P*ρ|ωρ(*t*)|, *x*4*m*(*t*) = 0.01*P*ρ|ϕρ(*t*)|,

$$\mathbf{x}\_{\mathsf{Sm}}(t) = 0.01 P\_{\mathsf{P}} |\boldsymbol{\delta}\_{\boldsymbol{\mathcal{P}}}(t)|\ldots$$

Here, the variable *P*<sup>ρ</sup> determines the relative width of the box *Bx* compared to the current values of the reference signals. This variable was increased to such a value that the condition (44) was met. The obtained value *P*<sup>ρ</sup> = 25% seems to be still admissible in the range of (43), (44) [4,5], determining the region of the local uniform asymptotic stability for the reference motion, which is realized by tracking controllers (35), (41).

#### **4. Discussion**

The main goal of this work was to propose constructive methods for marine tracking controllers' design taking into account the real conditions of a vessel's motions. We focused our main attention on a situation where the rudders' deflections and the forward speed are presented by initially given reference signals to be realized using tracking feedback controls.

This problem can be solved using different popular optimization approaches (Bellman's theory, MPC technique, sliding mode control, etc.). Nevertheless, we propose a new specific method for tracking controllers' design, which ensures the desirable reference motion of the vessel along the forward speed and heading angle.

This method is based on the optimal damping concept, which has certain advantages related to the practical requirements for the dynamic features of a closed-loop connection. The main advantage of the aforementioned approach is that the numerical solution of the optimization problems is essentially simplified. In contrast to well-known approaches [1–4], we applied OD tracking feedback as a control law with special features that allow it to be adjusted and implemented in real-time regime of motion.

The main result of this study is the development of the optimal damping concept to ensure its practical applicability and effectiveness that is illustrated by a controller design for a transport ship.

The investigations presented above could be further developed to consider the robust features of the tracking control laws, information about the measurement noise, and the presence of transport delays [24]. The results of the executed research could also be implemented to provide desirable reference motion of the dynamic positioning systems for marine vessels [26,27]. Certain attention may be given to the multipurpose control laws applications [22–25]. As for direct development of this study, it is possible to also use OD controller for the vertical rudders. The extension of the proposed approach to various robotic systems is also of considerable interest. The scope of the proposed approach may additionally include remotely operated vehicles [11] and offshore structures [12].

**Funding:** This research was funded by the Russian Foundation for Basic Research (RFBR) controlled by the Government of Russian Federation, research project number 20-07-00531.

**Acknowledgments:** The author expresses his sincere gratitude to A.P. Zhabko for the interest to the topic of this study and for the high-level professional suggestions.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **References**

