船舶运动具有大时滞、大惯性、强非线性等特点，同时由于航行条件和环境干扰的不确定性导致船舶运动模型具有严重的不确定性。因此，船舶航向控制成为非线性系统控制领域中的难点和热点问题^[1-2]，近年来不确定船舶航向系统控制问题取得了丰富的研究成果^[3-5]。现有大多数研究仅考虑了系统的稳定性，而良好的控制系统还应考虑系统的节能问题。神经网络和模糊逻辑系统作为通用逼近器，可以用来很好地辨别非线性系统的不确定性，学者将这2种方法应用到不确定船舶航向非系统的控制设计中，文献[6-7]采用两个神经网络，执行网络和评价网络逼近代价函数和控制律，分别基于近似动态规划算法和启发式动态规划算法^[8]，解决了船舶航向最优控制问题，但对于期望航向跟踪时间过长，并缺乏对船舶航向实际性能需求的考虑。本文研究了非线性离散时间船舶航向系统控制问题，利用航向跟踪误差设计效用函数，并定义策略效用函数。根据贝尔曼最优控制原理得到贝尔曼方程，设计代价函数。采用梯度下降法最小化代价函数，得到最优评价信号。将最优评价信号应用于控制器设计，不仅能够加快系统追踪速度，而且能够实现最优控制。

1 船舶航向运动数学模型

考虑船舶稳态回转非线性特性，本文采用如下船舶航向非线性系统数学模型：

$ \ddot \phi + \frac{1}{T}H(\dot \phi ) = \frac{K}{T}\delta $

(1)

式中：ϕ为航向角；δ为舵角；K是船舶转回性指数；T为船舶追随性指数；H($\dot \phi $)为$\dot \phi $的非线性函数可近似表示为$H\left( {\dot \phi } \right) = {a_1}\dot \phi + {a_2}{\dot \phi ^3} + {a_3}{\dot \phi ^5} + \ldots $，其中a_i(i=1, 2, 3, …)为常数。定义状态变量x₁=ϕ，x₂=${\dot \phi }$，u=δ，则式(1)经离散化，得到船舶航向非线性离散系统数学模型：

$ \left\{ {\begin{array}{*{20}{l}} {{x_1}(k + 1) = {x_2}(k)}\\ {{x_2}(k + 1) = {f_2}\left( {{x_2}(k)} \right) + pu(k)}\\ {{y_k} = {x_1}(k)} \end{array}} \right. $

(2)

式中：x_i(i=1, 2)为系统的状态；u(k)和y_k分别为系统的输入和输出；p=K/T为系统控制增益；函数f₂(x₂(k))=(-1/T)H(x₂(k))。

本文的控制目标是设计一个自适应模糊最优控制器u(k)，使式(2)的输出信号y_k能够跟踪参考信号y_d(k)，闭环系统半全局一致最终有界的，且跟踪误差收敛到零的一个小邻域内。控制系统共有2个假设条件：1)系统状态变量x₁(k), x₂(k)可测得；2)系统的跟踪参考信号为y_d(k)∈Ω_y, k>0，是已知并且有界的。

2 模糊逻辑系统及评价信号 2.1 模糊逻辑系统

本文采用具有线性可调参数的第一类模糊逻辑系统逼近未知连续函数^[9]，F_li和G_l为对应于模糊隶属函数μ_Fi^l(x_i)和的模糊集合，N是模糊规则数。通过模糊化和解模糊化方法，模糊逻辑系统为：

$ y(x) = \frac{{\sum\limits_{l = 1}^N {{{\bar y}_l}} \prod\limits_{i = 1}^n {{\mu _{F_i^l}}} \left( {{x_i}} \right)}}{{\sum\limits_{l = 1}^N {\left[ {\prod\limits_{i = 1}^n {{\mu _{F_i^l}}} \left( {{x_i}} \right)} \right]} }} $

(3)

其中${\bar y_l} = \mathop {{\rm{max}}}\limits_{y \in {\mathit{\pmb{R}}}} {\mu _{{G_l}}}\left( y \right)$。

定义模糊基函数为：

$ {\zeta _l} = \frac{{\prod\limits_{i = 1}^n {{\mu _{{F_{{l_i}}}}}} \left( {{x_i}} \right)}}{{\sum\limits_{l = 1}^N {\left( {\prod\limits_{i = 1}^n {{\mu _{{F_{{l_i}}}}}} \left( {{x_i}} \right)} \right)} }} $

(4)

由此，令${\mathit{\pmb{\theta}} ^{\rm{T}}} = \left[ {{{\bar y}_1},{{\bar y}_2}, \ldots ,{{\bar y}_N}} \right]{\rm{ = }}\left[ {{\theta _1},{\theta _2}, \ldots ,{\theta _N}} \right]$和$\zeta \left( x \right) = {\left[ {{\zeta _1}\left( x \right), {\zeta _2}\left( x \right), \ldots , {\zeta _N}\left( x \right)} \right]^{\rm{T}}}$，则模糊逻辑系统(3)可进一步表示为：

$ y(x) = {\mathit{\boldsymbol{\theta }}^{\rm{T}}}\zeta (x) $

(5)

引理1^[10]   f(x)是定义在闭集Ω内的连续函数，对任意给定的常数ε>0，则存在一个模糊逻辑系统(5)，使得:

$ \mathop {\sup }\limits_{x \in U} \left| {f(x) - {\mathit{\boldsymbol{\theta }}^{\rm{T}}}\zeta (x)} \right| \le \varepsilon $

(6)

基于引理1，函数f(x)可表示为^[11]：

$ f(x) = {\mathit{\boldsymbol{\theta }}^{*{\rm{T}}}}\zeta (x) + \mathit{\boldsymbol{\varepsilon }}(x) $

(7)

式中：θ^*和ε(x)代表最优参数向量和最小模糊逼近误差，满足$\left\| {{{\mathit{\pmb{\theta}}} ^*}} \right\| \le \bar \theta $和$\left\| {{\mathit{\pmb{\varepsilon}}} \left( x \right)} \right\| \le \bar \varepsilon $，其中常数$\bar \theta > 0$，$\bar \varepsilon > 0$。

2.2 最优评价信号

令追踪误差为e₁(k)=x₁(k)-y_d(k)，滤波误差:

$ e(k) = \mathit{\boldsymbol{ \boldsymbol{\varLambda} \bar e}}(k) $

(8)

式中$\bar {\mathit{\pmb{e}}}\left( k \right) = {[{e_1}\left( k \right), {e_1}\left( {k - 1} \right), \ldots , {e_1}\left( {k - n + 1} \right)]^{\rm{T}}}$，$\mathit{\pmb{\Lambda}} = [1, {\lambda _1}, \ldots , {\lambda _{n - 1}}]$，λ₁, …, λ_n为常数^[11]。通过模糊逻辑系统逼近未知非线性函数，利用贝尔曼方程设计目标函数并得到最优评价信号。基于滤波误差e(k)，效用函数q(k)为^[12]：

$ q(k) = \left\{ {\begin{array}{*{20}{l}} {0,}&{\left\| {\mathit{\boldsymbol{e}}(k)} \right\| \le c}\\ {1,}&{\left\| {\mathit{\boldsymbol{e}}(k)} \right\| > c} \end{array}} \right. $

(9)

式中c∈R为阈值，q(k)=0表示跟踪效果良好，q(k)=1代表跟踪效果不佳。

令策略效用函数C(k)为：

$ C(k) = q(k + 1) + \beta q(k + 2) + {\beta ^2}q(k + 3) + \cdots $

(10)

式中β>0是权重系数。

利用模糊逻辑系统逼近策略效用函数，则

$ C(k) = \mathit{\boldsymbol{\theta }}_c^{*{\rm{T}}}{\zeta _c}(k) + \mathit{\boldsymbol{\varepsilon }}_c^*(k) $

(11)

其中θ_c^*为最优模糊参数，${\hat {\mathit{\pmb{\theta}}} _c}\left( k \right)$为θ_c^*的估计值。

根据文献[12]和[15]，设贝尔曼方程

$ {e_c}(k) = \beta \hat C(k) - [\hat C(k - 1) - q(k)] $

(12)

式中：$\hat C(k) = \hat {\mathit{\pmb{\theta}}} _c^{\rm{T}}(k){\zeta _c}(k)$是C(k)的估计。

目标函数为E_c(k)=e_c²(k)/2，基于梯度下降规则，${\hat \theta _c}(k)$的自适应律为：

$ {\hat \theta _c}(k + 1) = {\hat \theta _c}(k) + \Delta {\hat \theta _c}(k) $

(13)

其中，$\Delta {{\hat \theta }_c}(k) = {\gamma _c}\left[ { - \frac{{\partial {E_c}(k)}}{{\partial {{\hat \theta }_c}(k)}}} \right] = {\gamma _c}\left[ { - \frac{{\partial {E_c}(k)}}{{\partial {e_c}(k)}}\frac{{\partial {e_c}(k)}}{{\partial {{\hat C}_c}(k)}}\frac{{\partial {{\hat C}_c}(k)}}{{\partial {{\hat \theta }_c}(k)}}} \right]$，学习参数γ_c>0。将式(12)代入式(13)，则

$ \begin{array}{*{35}{l}} {{{\hat{\theta }}}_{c}}(k+1)={{{\hat{\theta }}}_{c}}(k)-{{\gamma }_{c}}B{{\zeta }_{c}}(k)\cdot \\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ [\beta \hat{C}(k)-(\hat{C}(k-1)-q(k))] \\ \end{array} $

(14)

式(14)两边分别减去θ_c^*可得

$ \begin{array}{l} {{\tilde \theta }_c}(k + 1) = {{\tilde \theta }_c}(k) - {\gamma _c}\beta {\zeta _c}(k) \cdot \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;[\beta \hat C(k) - (\tilde C(k - 1) - q(k))] \end{array} $

(15)

式中${{\tilde \theta }_c}\left( k \right) = {{\hat \theta }_c}\left( k \right) - \theta _c^*$，令${\omega _c}\left( k \right) = {{\tilde \theta }_c}^T\left( k \right){\zeta _c}\left( k \right)$。

3 控制器设计和稳定性分析

船舶航向系统作为典型的严格反馈形式的非线性离散系统，基于自适应控制方法的控制器设计过程中会产生非因果关系问题，进而需要对式(2)等价变化^[13]：

$ \left\{ {\begin{array}{*{20}{l}} {{x_1}(k + 2) = {x_2}(k + 1)}\\ {{x_2}(k + 1) = {F_2}(k) + pu(k)}\\ {{y_k} = {x_1}(k)} \end{array}} \right. $

(16)

式中：F₂(k)=f₂(x₂(k))，针对转化后的船舶航向非线性离散系统，控制器设计及系统稳定性分析过程如下：

1) 设计虚拟控制函数。

根据系统式(16)，则:

$ {e_1}(k + 2) = {x_1}(k + 2) - {y_d}(k + 2) $

(17)

引入误差变量e₂(k+1)=x₂(k+1)-α₁(k)，设计虚拟控制函数α₁(k)为：

$ {\alpha _1}(k) = {y_d}(k + 2) $

(18)

由此可得：

$ {e_1}(k + 2) = {e_2}(k + 1) $

(19)

2) 设计优化控制器。

由e₂(k+1)=x₂(k+1)-α₁(k)，则：

$ {e_2}(k + 1) = {F_2}(k) + pu(k) - {\alpha _1}(k) $

(20)

定义函数

$ {U^*}(k) = - \left[ {{F_2}(k) - {\alpha _1}(k)} \right]/p $

(21)

利用模糊逻辑系统逼近效用函数U^*(k)为：

$ {U^*}(k) = {\mathit{\boldsymbol{\theta }}^{*{\rm{T}}}}\zeta (z(k)) + {\mathit{\boldsymbol{\varepsilon }}^*}(\mathit{\boldsymbol{z}}(k)) $

(22)

式中：z(k)=[x₁(k), x₂(k), y_d(k+2)]^T∈Ω_z，由式(21)~(20)可改写为：

$ {e_2}(k + 1) = p\left[ {u(k) - {\mathit{\boldsymbol{\theta }}^{*{\rm{T}}}}\zeta (\mathit{\boldsymbol{z}}(k)) - {\mathit{\boldsymbol{\varepsilon }}^*}(\mathit{\boldsymbol{z}}(k))} \right] $

(23)

实际控制输入u(k)为：

$ u(k) = {\mathit{\boldsymbol{\widehat \theta }}^{\rm{T}}}(k)\zeta (\mathit{\boldsymbol{z}}(k)) $

(24)

基于式(24)，(23)可改写为

$ \begin{array}{l} {e_2}(k + 1) = p\left[ {{{\mathit{\boldsymbol{\widehat \theta }}}^{\rm{T}}}(k)\zeta (z(k)) - } \right.\\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\left. {{\mathit{\boldsymbol{\theta }}^{*{\rm{T}}}}(k)\zeta (\mathit{\boldsymbol{z}}(k)) - {{\bf{ \pmb{\mathsf{ ε}} }}^*}(\mathit{\boldsymbol{z}}(k))} \right]\\ \;\;\;\;\;\;\;\;\;\;\;\;\;\; = p\left[ {\omega (k) - {{\bf{ \pmb{\mathsf{ ε}} }}^*}(\mathit{\boldsymbol{z}}(k))} \right] \end{array} $

(25)

式中$\omega \left( k \right) = {\tilde {\mathit{\pmb{\theta}}} ^{\rm{T}}}\left( k \right){\rm{ }}\zeta \left( {z\left( k \right)} \right)$，$\tilde \theta \left( k \right) = \hat \theta \left( k \right) - {\theta ^*}$。

用C_d(k)表示期望的策略效用函数，设：

$ {e_a}(k) = {\mathit{\boldsymbol{\widehat \theta }}^{\rm{T}}}(k)\zeta (\mathit{\boldsymbol{z}}(k)) + \left[ {\hat C(k) - {C_d}(k)} \right] $

(26)

通常，C_d(k)期望值为零代表系统输出具有较好的跟踪性能^[16]，故式(26)可表示为：

$ {e_a}(k) = {\mathit{\boldsymbol{\widehat \theta }}^{\rm{T}}}(k)\zeta (\mathit{\boldsymbol{z}}(k)) + \hat C(k) $

(27)

设目标函数E_a(k)=e_a²(k)/2，$\hat \theta \left( k \right)$的自适应律为：

$ \mathit{\boldsymbol{\widehat \theta }}(k + 1) = \mathit{\boldsymbol{\widehat \theta }}(k) + \Delta \mathit{\boldsymbol{\widehat \theta }}(k) $

(28)

式中$\Delta \hat \theta \left( k \right) = \gamma \left[ { - \frac{{\partial {E_a}\left( k \right)}}{{\partial \hat \theta \left( k \right)}}} \right]$，学习参数γ>0，则

$ \begin{array}{l} \hat \theta (k + 1) = \hat \theta (k) - \gamma \zeta (\mathit{\boldsymbol{z}}(k)) \cdot \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\left[ {{{\mathit{\boldsymbol{\widehat \theta }}}^{\rm{T}}}(k)\zeta (\mathit{\boldsymbol{z}}(k)) + \dot C(k)} \right] \end{array} $

(29)

定理1：  对于船舶航向离散非线性系统式(16)，实际控制器式(24)，自适应律式(29)构成的闭环系统所有信号半全局一致最终有界，通过调整设计参数使得追踪误差收敛到零的一个小邻域内。

证明：  首先，选取闭环系统李雅普诺夫预选函数V(k)的形式为：

$ V(k) = {V_1}(k) + {V_2}(k) + {V_3}(k) + {V_4}(k) + {V_5}(k) $

(30)

其中，V₁(k)=τ_e1e₁²(k)/2，V₂(k)=τ_e2e₂²(k)/2，${V_3}\left( k \right) = {\rm{ }}{\tau _\theta }{{\tilde {\mathit{\pmb{\theta}}} }^{\rm T}}\left( k \right)\tilde \theta \left( k \right)/\gamma $，${V_4}\left( k \right) = {\tau _c}\tilde {\mathit{\pmb{\theta}}} _c^{\rm T}\left( k \right)_c^{\rm T}{{\tilde {\mathit{\pmb{\theta}}} }_c}\left( k \right)/{\gamma _c}$，V₅(k)=2τ_cω_c²(k-1)。

式中：τ_e1、τ_e2、τ_θ、τ_c为正常数。V₁(k)的差分为：

$ \begin{array}{l} \Delta {V_1}(k) = {\tau _{{e_1}}}\left[ {e_1^2(k + 1) - e_1^2(k)} \right]/2 = \\ \;\;\;\;\;\;\;\;\;\;\;\;\;{\tau _{{e_1}}}\left[ {e_2^2(k) - e_1^2(k)} \right]/2 \end{array} $

(31)

V₂(k)的差分为：

$ \begin{array}{l} \Delta {V_2}(k) = {\tau _{{e_2}}}\left[ {e_2^2(k + 1) - e_2^2(k)} \right]/2 = \\ \;\;\;\;\;\;\;\;\;\;\;\;\;{\tau _{{e_2}}}\left[ {p{{\left( {\omega (k) - {\varepsilon ^*}(z(k))} \right]}^2}/2} \right. - \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;{\tau _{{e_2}}}e_2^2(k)/2 \end{array} $

(32)

基于Cauchy-Schwarz不等式：

$ {\left( {{a_1} + \cdots + {a_n}} \right)^2} \le n\left( {a_1^2 + \cdots + a_n^2} \right) $

(33)

式(32)变为:

$ \Delta {V_2}(k) \le {\tau _{{e_2}}}{p^2}\left( {{\omega ^2}(k) + {{\bar \varepsilon }^2}} \right) - {\tau _{{e_2}}}e_2^2(k)/2 $

(34)

V₃(k)的差分为：

$ \begin{array}{*{20}{c}} {\Delta {V_3}(k) = {\tau _\theta }{{\tilde \theta }^{\rm{T}}}(k + 1)\tilde \theta (k + 1)/\gamma - }\\ {{\tau _\theta }{{\tilde \theta }^{\rm{T}}}(k)\tilde \theta (k)/\gamma } \end{array} $

(35)

式(29)两边减去θ^*可得：

$ \begin{array}{l} \tilde \theta (k + 1) = \tilde \theta (k) - \gamma \zeta (z(k)) \cdot \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\left[ {{{\mathit{\boldsymbol{\widehat \theta }}}^{\rm{T}}}(k)\zeta (\mathit{\boldsymbol{z}}(k)) + \hat C(k)} \right] \end{array} $

(36)

由(36)，则式(35)变为:

$ \begin{array}{l} \Delta {V_3}(k) = {\tau _\theta }[\tilde \theta (k) - \gamma \zeta (z(k)) \cdot \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;{\left. {\left( {{{\mathit{\boldsymbol{\widehat \theta }}}^{\rm{T}}}(k)\zeta (\mathit{\boldsymbol{z}}(k)) + \hat C(k)} \right)} \right]^{\rm{T}}} \cdot \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;[\mathit{\widetilde \theta }(k) - \gamma \zeta (z(k)) \cdot \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\left. {\left( {{{\mathit{\boldsymbol{\widehat \theta }}}^{\rm{T}}}(k)\zeta (\mathit{\boldsymbol{z}}(k)) + \hat C(k)} \right)} \right]/\gamma - \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;{\tau _\theta }{\mathit{\boldsymbol{\widetilde \theta }}^{\rm{T}}}(k)\mathit{\boldsymbol{\widetilde \theta }}(k)/\gamma = \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;{\tau _\theta }\gamma {\left\| {\zeta \left( {z\left( k \right)} \right)} \right\|^2}\left[ {{{\mathit{\boldsymbol{\widehat \theta }}}^{\rm{T}}}(k)} \right. \cdot \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\zeta (z(k)) + \hat C(k){]^2} - 2{\tau _\theta }\mathit{\boldsymbol{\omega }}(k) \cdot \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\left[ {{{\mathit{\boldsymbol{\widehat \theta }}}^{\rm{T}}}(k)\zeta (\mathit{\boldsymbol{z}}(k)) + \hat C(k)} \right] = \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\begin{array}{*{20}{l}} { - {\tau _\theta }\left[ {1 - \gamma {{\left\| {\zeta (\mathit{\boldsymbol{z}}(k))} \right\|}^2}} \right]}\\ {{{\left[ {{{\mathit{\boldsymbol{\widehat \theta }}}^{\rm{T}}}(k)\zeta (\mathit{\boldsymbol{z}}(k)) + \hat C(k)} \right]}^2} + }\\ {{\tau _\theta }\left[ {{\mathit{\boldsymbol{\theta }}^{*{\rm{T}}}}\zeta (\mathit{\boldsymbol{z}}(k)) + \hat C(k)} \right] \cdot }\\ {\left[ {{{\mathit{\boldsymbol{\widehat \theta }}}^{\rm{T}}}(k)\zeta (\mathit{\boldsymbol{z}}(k)) + \hat C(k)} \right] - }\\ {{\tau _\theta }\omega (k)\left[ {{{\mathit{\boldsymbol{\widehat \theta }}}^{\rm{T}}}(k)\zeta (\mathit{\boldsymbol{z}}(k)) + \hat C(k)} \right]} \end{array} \end{array} $

(37)

因为‖ζ(z(k))‖²≤1，则:

$ \begin{array}{l} \Delta {V_3}(k) \le - {\tau _\theta }(1 - \gamma )\left[ {{{\mathit{\boldsymbol{\widehat \theta }}}^{\rm{T}}}(k)\zeta (\mathit{\boldsymbol{z}}(k)) + } \right.\\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\hat C(k){]^2} + {\tau _\theta }\left[ {{\mathit{\boldsymbol{\theta }}^{*{\rm{T}}}}\zeta (\mathit{\boldsymbol{z}}(k)) + } \right.\\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\hat C(k){]^2} - {\tau _\theta }{\mathit{\boldsymbol{\omega }}^2}(k) \end{array} $

(38)

利用式(33)、式(38)变为：

$ \begin{array}{l} \Delta {V_3}(k) \le - {\tau _\theta }(1 - \gamma )\left[ {{{\mathit{\boldsymbol{\widehat \theta }}}^{\rm{T}}}(k)\zeta (\mathit{\boldsymbol{z}}(k)) + } \right.\\ \;\;\;\;\;\;\;\;\;\;\;\;\;\hat C(k){]^2} + {\tau _\theta }\left[ {{\mathit{\boldsymbol{\theta }}^{*{\rm{T}}}}\zeta (\mathit{\boldsymbol{z}}(k)) + } \right.\\ \;\;\;\;\;\;\;\;\;\;\;\;\;{\left. {\mathit{\boldsymbol{\widehat \theta }}_c^{\rm{T}}(k){\zeta _c}(k)} \right]^2} - {\tau _\theta }{\mathit{\boldsymbol{\omega }}^2}(k) \le \\ \;\;\;\;\;\;\;\;\;\;\;\;\; - {\tau _\theta }(1 - \gamma )\left[ {{{\hat \theta }^{\rm{T}}}(k)\zeta (\mathit{\boldsymbol{z}}(k))} \right. + \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\hat C(k){]^2} + 2{\tau _\theta }\omega _c^2(k) - \\ \;\;\;\;\;\;\;\;\;\;\;\;\;{\tau _\theta }{\omega ^2}(k) + 2{\tau _\theta }{\left( {\bar \theta + {{\bar \theta }_c}} \right)^2} \end{array} $

(39)

V₄(k)的差分为:

$ \begin{array}{l} \Delta {V_4}(k) = {\tau _c}\left[ {\mathit{\boldsymbol{\widetilde \theta }}_c^{\rm{T}}(k + 1){{\tilde \theta }_c}(k + 1) - } \right.\\ \;\;\;\;\;\;\;\;\;\;\;\left. {\mathit{\boldsymbol{\widetilde \theta }}_c^{\rm{T}}(k){{\tilde \theta }_c}(k)} \right]/{\gamma _c} = \\ \;\;\;\;\;\;\;\;\;\;\;{\tau _c}\left[ {{{\tilde \theta }_c}(k) - {\gamma _c}\beta {\zeta _c}(k)} \right. \cdot \\ \;\;\;\;\;\;\;\;\;\;\;(\beta \tilde C(k) + q(k) - \hat C(k - 1)){]^{\rm{T}}} \cdot \\ \;\;\;\;\;\;\;\;\;\;\;\left[ {{{\tilde \theta }_c}(k) - {\gamma _c}\beta {\zeta _c}(k)} \right. \cdot \\ \;\;\;\;\;\;\;\;\;\;\;\begin{array}{*{20}{l}} {(\beta \hat C(k) + q(k) - \hat C(k - 1))]/{\gamma _c} - }\\ {{\tau _c}\tilde \theta _c^T(k){{\tilde \theta }_c}(k)/{\gamma _c} = }\\ { - {\tau _c}\left( {1 - {\gamma _c}{\beta ^2}{{\left\| {{\zeta _c}(k)} \right\|}^2}} \right)[\beta \hat C(k) + }\\ {q(k) - \hat C(k - 1){]^2} - 2{\tau _c}\beta {\omega _c}(k)}\\ {[\beta \hat C(k) + q(k) - \hat C(k - 1)] + }\\ {{\tau _c}{{[\beta \hat C(k) + q(k) - \hat C(k - 1)]}^2}} \end{array} \end{array} $

(40)

因为‖ζ_c(k)‖²≤1，式(40)变为:

$ \begin{array}{l} \Delta {V_4}\left( k \right) \le - {\tau _c}\left( {1 - {\gamma _c}{\beta ^2}} \right)\left[ {\beta \hat C\left( k \right) + q\left( k \right) - } \right.\\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;{\left. {\hat C\left( {k - 1} \right)} \right]^2} + {\tau _c}\left[ {\beta \mathit{\boldsymbol{\theta }}_c^{ * {\rm{T}}}{\zeta _c}\left( k \right) + } \right.\\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;q\left( k \right) - {\omega _c}\left( {k - 1} \right) - \\ {\left. {\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\mathit{\boldsymbol{\theta }}_c^{*{\rm{T}}}{\zeta _c}(k - 1)} \right]^2} - {\tau _c}{\beta ^2}\omega _c^2(k) \end{array} $

(41)

利用式(33)、式(41)变为:

$ \begin{array}{l} \Delta {V_4}(k) \le - {\tau _c}\left( {1 - {\gamma _c}{\beta ^2}} \right)[\beta \hat C(k) + \\ \;\;\;\;\;\;\;\;\;\;\;\;\;q(k) - \hat C(k - 1){]^2} - \\ \;\;\;\;\;\;\;\;\;\;\;\;\;{\tau _c}{\beta ^2}\omega _c^2(k) + 2{\tau _c}\omega _c^2(k - 1) + \\ \;\;\;\;\;\;\;\;\;\;\;\;\;2{\tau _c}{\left( {_c^{\bar \theta }\left( {1 + \beta } \right) + 1} \right)^2} \end{array} $

(42)

V₅(k)的差分为：

$ \Delta {V_5}(k) = 2{\tau _c}\omega _c^2(k) - 2{\tau _c}\omega _c^2(k - 1) $

(43)

结合式(31)、式(34)、式(39)、式(42)和(43)则：

$ \begin{array}{l} \Delta V(k) \le - {\tau _{{e_1}}}e_1^2(k)/2 - \left( {{\tau _{{e_2}}} - {\tau _{{e_1}}}} \right)e_2^2(k)/2 - \\ \;\;\;\;\;\;\;\;\;\;\;\;\;{\tau _\theta }(1 - \gamma ){\left[ {{{\mathit{\boldsymbol{\widehat \theta }}}^{\rm{T}}}(k)\zeta (z(k)) + \hat C(k)} \right]^2} - \\ \;\;\;\;\;\;\;\;\;\;\;\;\;{\tau _c}\left( {1 - {\gamma _c}{\beta ^2}} \right)[\beta \hat C(k) + q(k) - \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\begin{array}{*{20}{l}} {\hat C(k - 1){]^2} - \left( {{\tau _\theta } - {\tau _{{e_2}}}{p^2}} \right){\omega ^2}(k) - }\\ {\left( {{\tau _c}{\beta ^2} - 2{\tau _\theta } - 2{\tau _c}} \right)\omega _c^2(k) + D} \end{array} \end{array} $

(44)

式中$D = 2{\tau _\theta }{(\bar \theta + {\bar \theta _c})^2} + {\tau _{{e_2}}}{p^2}{\bar \varepsilon ^2} + 2{\tau _c}{[\left( {1 + \beta } \right){\bar \theta _c} + 1]^2}$。

选取0 < β² < 1/γ_c，0 < γ < 1则，

$ \begin{array}{l} \Delta V(k) \le - {\tau _{{e_1}}}e_1^2(k)/2 - \left( {{\tau _{{e_2}}} - {\tau _{{e_1}}}} \right)e_2^2(k)/2 - \\ \;\;\;\;\;\;\;\;\;\;\;\;\left( {{\tau _\theta } - {\tau _{{e_2}}}{p^2}} \right){\omega ^2}(k) - \\ \;\;\;\;\;\;\;\;\;\;\;\;\left( {{\tau _c}{\beta ^2} - 2{\tau _\theta } - 2{\tau _c}} \right)\omega _c^2(k) + D \end{array} $

(45)

设计参数满足τ_e1>0，τ_e2>τ_e1，τ_θ>τ_e2p²，τ_cβ²>2τ_θ+2τ_c。类似文献^[12-16]的分析，参数满足$\left| {{e_1}\left( k \right)} \right| > \sqrt D /\sqrt {{\tau _{e1}}/2} $，或$\left| {{e_2}\left( k \right)} \right| > \sqrt D /\sqrt {{\tau _{{e_2}}}/2 - {\tau _{{e_1}}}/2} $，或$\left| {\omega \left(k \right)} \right| > \sqrt D /\sqrt {{\tau _\theta } - {\tau _{{e_2}}}{p^2}} $，或$\left| {{\omega _c}\left( k \right)} \right| > \sqrt D /\sqrt {{\tau _c}{\beta ^2} - 2{\tau _\theta } - 2{\tau _c}} $任意一个不等式条件之一，则ΔV(k) < 0，根据现有Lyapunov理论可以证明，闭环系统内所有信号是半全局一致最终有界的，系统跟踪误差收敛到以零为中心的邻域内。

4 有效性实例仿真研究

利用Matlab进行计算机仿真，以大连海事大学远洋实习船“育龙”轮为例，通过与文献[14]算法的对比研究，验证本文控制算法的有效性。跟踪信号选取能够代表实际性能要求的数学模型：

$ \begin{array}{*{20}{l}} {{\phi _m}(k + 2) + 0.1{\phi _m}(k + 1) + 0.0025{\phi _m}(k) = }\\ {\;\;\;0.0025{\phi _r}(k)} \end{array} $

(46)

式中：ϕ_m代表船舶航向的期望系统性能；ϕ_r(k)=(sign(sin(πk/500))+1)π/12是一个经过处理的输入信号，其取值为0°~30°，周期为500 s。

已知“育龙”轮船长126.0 m.，船宽20.8 m, 满载吃水8.0 m, 方形系数为0.681，船速7.7 m/s。通过计算得船舶航向离散非线性系统数学模型参数a₁=1, a₂=30, K=0.478, T=216。

在仿真中，船舶航向系统的初始条件选取为：γ=0.05的初始状态为5°(即航向偏差为+5°)，γ=0.05的初始状态为0，待设计参数选取γ=0.05，γ_c=0.01，β=0.05。模糊隶属函数选取为：

$ \begin{array}{l} \begin{array}{*{20}{l}} {{\mu _{{F_{i1}}}}(x) = \exp \left[ { - {{(x + 2)}^2} + 3} \right]}\\ {{\mu _{{F_{i2}}}}(x) = \exp \left[ { - {{(x + 1)}^2} + 3} \right]} \end{array}\\ \begin{array}{*{20}{l}} {{\mu _{{F_{i3}}}}(x) = \exp \left[ { - {{(x)}^2} + 3} \right]}\\ {{\mu _{{F_{i4}}}}(x) = \exp \left[ { - {{(x - 1)}^2} + 3} \right]} \end{array}\\ {\mu _{{F_{i5}}}}(x) = \exp \left[ { - {{(x - 2)}^2} + 3} \right] \end{array} $

(47)

其中x=[x₁(k), x₂(k), y_d(k+2)]^T，考虑风浪流等外界干扰的影响，本文在仿真中加了干扰项0.1sint。利用Matlab进行计算机仿真研究，图 1为船舶航向非线性离散系统自适应模糊最优控制算法仿真图。

由图 1可见，本文相比文献[14]，在考虑性系统节能的前提下，对于期望航向跟踪时间较短，兼顾了船舶航向实际性能需求；在跟踪过程中，误差趋于以零为中心的一个小邻域内，表明跟踪效果良好；控制舵角的历时曲线无较大超调量并趋于稳定，表明满足要求；自适应参数的历时曲线显示其为有界的。以上仿真结果分析验证了本文所提算法的有效性。

5 结论

1) 与现有船舶航向离散非线性系统的控制方法相比较，本文考虑控制器节能问题，较好的保持了航向控制系统性能，但同时也使系统的稳定性分析更为复杂和困难。

2) 经研究得出，算法在实现最优控制的同时，确保了闭环系统达到半全局最终一致稳定，也使系统的跟踪误差收敛到以零为中心的邻域内，是解决考虑节能问题的船舶航向控制问题的一个更好的方法。

工程应用中较难达到系统状态信息完全已知的情况，本文没有解决系统状态不可测的问题，下一步研究将从容错控制和状态约束等方面，基于状态观测器，探讨船舶航向控制输出反馈控制问题。