无人船(umanned surface vessels，USV)是船舶未来发展的趋势，近年，随着物联网、大数据、云计算、人工智能等新技术的突飞猛进，船舶自动化的水平不断提高，无人船舶的实现有了科技的支撑，无人驾驶船舶航行于全球的实现有了可能性^[1]。无人船舶作为一个软硬件相结合、功能强大的复杂系统，安全且可靠的航行需要船载硬件容错管理、传感器信息融合以及规划控制等多个功能模块的协同配合工作。规划控制作为一种可实现无人船舶自主航行的系统化控制策略，在广义上可由路径规划、行为决策、动作规划以及反馈控制4个部分组成。其中航向反馈控制问题是实现无人船舶自主航行的最底层问题，也是规划控制策略中最基本的反馈控制执行问题，可直接决定规划控制策略的整体控制效果。由于船舶航向控制问题作为国内外船舶控制领域的重点研究课题，有非常好的研究基础进行借鉴。例如双滤波器的航向保持控制方法^[2]、模糊神经网络船舶航向自适应控制方法^[3]、基于反步法设计的神经网络航向保持控制方法^[4]、基于最小二乘支持向量机的航向保持模糊控制方法^[5]。基于李雅普诺夫稳定性理论的自适应鲁棒航向保持控制方法^[6]，以及基于非线性反馈鲁棒控制策略^[7]、专家系统^[8]、滑模变结构控制^[9]的航向保持控制方法。上述研究仅局限于单纯船舶航向控制算法的研究，而且计算复杂度高的算法引入会导致航向控制算法在规划控制策略上层计算的复杂度爆炸性上升，甚至发生不同规划控制模块之间的冲突。因此设计一种计算复杂度低，结构简捷，控制效果和鲁棒性俱佳的航向反馈控制器，是实现无人船舶规划控制策略的关键问题。闭环增益成形算法^[10-12]结构简单、物理意义明确，非常适合作为规划控制策略中船舶的航向反馈控制模块的控制方案，因此本文通过理论分析和仿真试验，对基于闭环增益成形算法所设计的控制器最终应用于无人船的可行性进行了研究。

1 船舶运动数学模型

本文以响应型船舶运动模型为例来表述无人船的运动。考虑外界影响下各种改变舵角的操纵运动作为输出操纵运动对输入舵角的响应关系，并沿用Abkowitz的研究方案^[10]，把流体动力X、Y、N展开成Taylor级数时只保留一阶小量，并假设船舶处在等速直线运动这一平衡状态，只考虑船舶航向控制与横漂及转艏运动之间的强耦合性，得到如下方程:

$ \begin{array}{*{20}{c}} {\left[ {\begin{array}{*{20}{c}} {\left( {m - {Y_{\dot v}}} \right)}&{\left( {m{x_C} - {Y_{\dot r}}} \right)}\\ {\left( {m{x_C} - {N_{\dot v}}} \right)}&{\left( {{I_{zz}} - {N_{\dot r}}} \right)} \end{array}} \right]\left[ {\begin{array}{*{20}{c}} {\dot v}\\ {\dot r} \end{array}} \right] = }\\ {\left[ {\begin{array}{*{20}{c}} {{Y_v}}&{\left( {{Y_r} - m{u_0}} \right)}\\ {{N_v}}&{\left( {{N_r} - m{x_C}{u_0}} \right)} \end{array}} \right]\left[ {\begin{array}{*{20}{l}} v\\ r \end{array}} \right] + \left[ {\begin{array}{*{20}{c}} {{Y_\delta }}\\ {{N_\delta }} \end{array}} \right]\delta } \end{array} $

(1)

式中:v和r分别代表船舶的横漂速度和转艏角速度; δ为船舶的舵角; m为船舶质量; x_C为船舶质心到船舶运动坐标系原点的距离; I_zz为船舶转动惯性矩; Y和N分别表示船舶在横漂方向和绕质心转动的外力。如果将式(1)所示矩阵方程的第1行两端除以ρL³/2，第2行除以ρL⁴/2，并转化为无量纲流体动力导数，可得到如下标准的状态空间形式:

$ {{\mathit{\boldsymbol{\dot X}}}_{(3)}} = {\mathit{\boldsymbol{A}}_{(3)}}{\mathit{\boldsymbol{X}}_{(3)}} + {\mathit{\boldsymbol{B}}_{(3)}}\delta $

(2)

式中:系统状态向量X₍₃₎=[v r ψ]^T，其中$\dot{\psi}=r$，${\mathit{\boldsymbol{A}}_{(3)}} = \left[ {\begin{array}{*{20}{c}} {{a_{11}}}&{{a_{12}}}&0\\ {{a_{21}}}&{{a_{22}}}&0\\ 0&1&0 \end{array}} \right], {\mathit{\boldsymbol{B}}_{(3)}} = \left[ {\begin{array}{*{20}{c}} {{b_{11}}}\\ {{b_{21}}}\\ 0 \end{array}} \right]$，其中a₁₁、a₁₂、a₂₁、a₂₂、b₁₁、b₂₁一般可根据船舶的船型参数估算得出^[10]。为了完整地描述船舶航向运动，在式(2)中加入输出方程:

$ \psi = \mathit{\boldsymbol{C}}{\mathit{\boldsymbol{X}}_{(3)}} $

式中C=[0 0 1]，可将状态空间模型转化为传递函数形式为:

$ {G_{\psi \delta }}(s) = \mathit{\boldsymbol{C}}{(s\mathit{\boldsymbol{I}} - \mathit{\boldsymbol{A}})^{ - 1}}\mathit{\boldsymbol{B}} = \frac{{{K_0}\left( {{T_3}s + 1} \right)}}{{s\left( {{T_1}s + 1} \right)\left( {{T_2}s + 1} \right)}} $

令1/T₁T₂=a₁₁a₂₂-a₁₂a₂₁，(T₁+T₂)/T₁T₂=-(a₁₁+a₂₂)，1/T₃=(b₁₁a₂₁-b₂₁a₁₁)/b₂₁，K₀T₃/T₁T₂=b₂₁，由此得到3个船舶时间常数T₁、T₂、T₃以及系统增益系数K₀。

取舵角δ为被控对象的输入，艏向角ψ为被控对象输出，描述操舵对艏摇响应的传递函数为:

$ {G_{\psi \delta }}(s) = \frac{{\psi (s)}}{{\delta (s)}} = \frac{{{K_0}\left( {1 + {T_3}s} \right)}}{{s\left( {1 + {T_1}s} \right)\left( {1 + {T_2}s} \right)}} $

对于船舶这种大惯性大滞后的运载工具，其动态特征在低频段尤为重要，上式中令s=jω→0，且再由泰勒级数展开忽略高阶小量导出如下简化模型:

$ {G_{\psi \delta }}(s) = \frac{{{K_0}}}{{s\left( {{T_0}s + 1} \right)}} $

(3)

式中:K₀代表船舶舵增益系数; T₀代表船舶追随性时间常数，T₀=T₁+T₂-T₃。式(3)也可表示为:

$ {T_0}\ddot \psi + \dot \psi = {K_0}\delta $

由于船舶在航行中还受到风浪流等作用，因此对上述模型的修正如下:

$ {T_0}\ddot \psi + \dot \psi = {K_0}\delta + w $

(4)

式中w为有界干扰项，即‖w‖_∞≤ρ。

2 基于闭环增益成形算法的无人船航向保持控制器设计

闭环增益成形算法作为一种工程简化的H_∞鲁棒控制算法，可针对如图 1所示的标准反馈控制系统，运用混合灵敏度控制算法的结果，利用4个参数(最大奇异值、带宽频率、关门频率及闭环频谱峰值)构造出补灵敏度函数，再运用补灵敏度函数T与灵敏度函数S=1/(GK+1)的相关性(T+S=I)间接构造出灵敏度函数S，最终得到控制器K，可以说该算法就是利用对整个闭环系统补灵敏度函数的成形，从而获得与混合灵敏度算法相同的控制性能与鲁棒性能。

为方便对图 1所示标准反馈控制系统的信号跟踪问题进行分析，做如下假设:

假设1:闭环系统的工作频谱为低通，即系统的主要工作在低频区域，高频区域作为干扰存在;

假设2:闭环系统的频谱最大奇异值设置为1，即系统的频谱峰值为0，从而保证系统稳定、无静差地跟踪参考信号r;

假设3:闭环系统频谱的关门斜率可取-20 dB/dec，用来对外界干扰的抑制。

补充说明:假设3中的闭环系统频率关门斜率的选择可以根据实际需要，选择斜率更大的-40 dB/dec和-60 dB/dec，用以增强系统的抗干扰性能，但随着关门频率的变大，设计出的控制器阶数会相应变高，使得控制器的结构和实现更复杂。

在假设1~3下，整个闭环控制系统的控制效果将直接由系统的带宽频率1/T₁所决定，而且整个控制器的设计过程简单，所调参数的物理意义明确，是一种简捷控制算法^[11]。

针对式(3)所示船舶航向运动数学模型，在假设1~3的条件下，采用闭环增益成形算法进行控制器设计，并将闭环系统的带宽频率设置为1/T₁，则此时航向控制系统补灵敏度函数为:

$ \frac{{{G_{\psi \delta }}(s)K(s)}}{{1 + {G_{\psi \delta }}(s)K(s)}} = \frac{1}{{{T_1}s + 1}} $

(5)

将船舶航向模型式(3)代入式(5)，最终得到船舶航向控制器:

$ K = \frac{1}{{{K_0}{T_1}}} + \frac{{{T_0}}}{{{K_0}{T_1}}}s $

(6)

式(6)所设计的航向控制器实质上就是一种基于负反馈机制的比例-微分(proportional differential，PD)控制器，其继承了常规PD控制器结构简单、使用方便的特点，同时本文所提方法也解决了常规PD控制器参数整定繁杂、物理意义不明确的问题。这使无人船在反馈控制层的设计上可以在不降低控制性能的前提下，解决了计算复杂度上升以及与无人船规划控制策略中与其他模块冲突的问题。

3 控制器稳定性分析与证明

闭环增益成形算法一直缺乏完善的控制器稳定性分析与证明，本文采用Lyapunov稳定性理论对基于闭环增益成形算法所设计的船舶航向反馈控制器进行稳定性分析。

首先定义控制系统的状态变量X=[x₁ x₂]^T，其中:

$ \left\{ {\begin{array}{*{20}{l}} {{\mathit{\boldsymbol{x}}_1} = e = {\psi _r} - \psi }\\ {{\mathit{\boldsymbol{x}}_2} = {{\dot x}_1} = - \dot \psi } \end{array}} \right. $

(7)

根据式(6)可得到船舶航向控制系统的舵角δ的表达式为:

$ \delta = Ke = \frac{1}{{{K_0}{T_1}}}e + \frac{{{T_0}}}{{{K_0}{T_1}}}\dot e $

(8)

把式(8)代入式(4)中，则:

$ {T_0}\ddot \psi + \dot \psi = {K_0}\left( {\frac{1}{{{K_0}{T_1}}}e + \frac{{{T_0}}}{{{K_0}{T_1}}}\dot e} \right) + w $

(9)

根据式(7)、(9)整理可得被控系统的状态方程:

$ \mathit{\boldsymbol{\dot X}} = \mathit{\boldsymbol{AX}} + \mathit{\boldsymbol{Bw}} $

(10)

其中:$\mathit{\boldsymbol{A}} = \left[ {\begin{array}{*{20}{c}} 0&1\\ { - \frac{1}{{{T_0}{T_1}}}}&{ - \frac{{{T_1} + {T_0}}}{{{T_0}{T_1}}}} \end{array}} \right], \mathit{\boldsymbol{B}} = \left[ {\begin{array}{*{20}{c}} 0\\ { - \frac{1}{{{T_0}}}} \end{array}} \right]$。

针对式(10)，定义系统的Lyapunov函数为:

$ \mathit{\boldsymbol{V}}(x) = {\mathit{\boldsymbol{X}}^{\rm{T}}}\mathit{\boldsymbol{PX}} $

(11)

式中P为正定实对称矩阵，且P=P^T。

在不考虑式(10)中外界干扰项w，则可得到:

$ \mathit{\boldsymbol{\dot V}}(\mathit{\boldsymbol{x}}) = {\mathit{\boldsymbol{X}}^{\rm{T}}}\mathit{\boldsymbol{PAX}} + {(\mathit{\boldsymbol{AX}})^{\rm{T}}}\mathit{\boldsymbol{PX}} = {\mathit{\boldsymbol{X}}^{\rm{T}}}\left( {\mathit{\boldsymbol{PA}} + {\mathit{\boldsymbol{A}}^{\rm{T}}}\mathit{\boldsymbol{P}}} \right){\mathit{\boldsymbol{X}}}。$

根据Lyapunov稳定性定理可知，对任意给定的正定实对称矩阵Q，欲使系统在原点处大范围渐进稳定，需满足以下方程:

$ {\mathit{\boldsymbol{A}}^{\rm{T}}}\mathit{\boldsymbol{P}} + \mathit{\boldsymbol{PA}} = - \mathit{\boldsymbol{Q}} $

(12)

定义正定实对称矩阵$\mathit{\boldsymbol{Q}} = \mathit{\boldsymbol{I}} = \left[ {\begin{array}{*{20}{l}} 1&0\\ 0&1 \end{array}} \right]$，并构造正定实对称矩阵P为:

$ \mathit{\boldsymbol{P}} = \left[ {\begin{array}{*{20}{l}} {{P_{11}}}&{{P_{12}}}\\ {{P_{12}}}&{{P_{22}}} \end{array}} \right] $

(13)

将式(10)、(13)代入式(12)可得到:

$ \begin{array}{*{20}{c}} {\left[ {\begin{array}{*{20}{c}} 0&1\\ { - \frac{1}{{{T_0}{T_1}}}}&{ - \frac{{{T_1} + {T_0}}}{{{T_0}{T_1}}}} \end{array}} \right]\left[ {\begin{array}{*{20}{c}} {{P_{11}}}&{{P_{12}}}\\ {{P_{12}}}&{{P_{22}}} \end{array}} \right] + }\\ {\left[ {\begin{array}{*{20}{c}} {{P_{11}}}&{{P_{12}}}\\ {{P_{12}}}&{{P_{22}}} \end{array}} \right]\left[ {\begin{array}{*{20}{c}} 0&1\\ { - \frac{1}{{{T_0}{T_1}}}}&{ - \frac{{{T_1} + {T_0}}}{{{T_0}{T_1}}}} \end{array}} \right] = - \left[ {\begin{array}{*{20}{c}} 1&0\\ 0&1 \end{array}} \right]} \end{array} $

(14)

对矩阵方程(14)展开并求解可得:

$ {P_{11}} = \left( {{{\left( {{T_1} + {T_0}} \right)}^2} + \left( {1 + {T_0}{T_1}} \right)} \right)/\left( {2\left( {{T_0} + {T_1}} \right)} \right), $

$ {P_{12}} = {T_0}{T_1}/2, $

$ {P_{22}} = \left( {{T_0}{T_1}\left( {1 + {T_0}{T_1}} \right)} \right)/\left( {2\left( {{T_0} + {T_1}} \right)} \right) $

为了保证矩阵P的正定性，需要满足不等式P₁₁>0以及P₁₁P₂₂-P₁₂²>0，最终求解可得:

$ \left\{ {\begin{array}{*{20}{l}} {{T_1} > \frac{{ - 3{T_0} + \sqrt {5T_0^2 - 4} }}{2}}\\ {{T_0}{T_1} > 0} \end{array}} \right. $

(15)

在满足${T_0} > \sqrt {0.8} $条件下，只要保证T₁>0即可满足控制器在原点处的平衡状态大范围渐进稳定。

为进一步论证本文所提出控制器对系统外界干扰的鲁棒性能，考虑式(10)中的外界干扰项w，则如式(11)所示的Lyapunov函数的导数为:

$ \dot V(x) = - x_1^2 - x_2^2 - {T_1}w{x_1} - kw{x_2} $

其中，k=T₁(1+T₀T₁)/(T₁+T₀)。

根据Youth′s不等式，可知:

$ \left\{ {\begin{array}{*{20}{l}} { - wT{x_1} \le T_1^2x_1^2 + \frac{{{w^2}}}{4}}\\ { - kw{x_2} \le {k^2}x_2^2 + \frac{{{w^2}}}{4}} \end{array}} \right. $

则:

$ \begin{array}{*{20}{c}} {\dot V(x) = - \left( {1 - T_1^2} \right)x_1^2 - \left( {1 - T_1^2} \right)x_2^2 + \frac{{{w^2}}}{2} \le }\\ { - \left( {1 - T_1^2} \right)x_1^2 - \left( {1 - {k^2}} \right)x_2^2 + \frac{{{\rho ^2}}}{2}} \end{array} $

定义a=max(T₁, k)时，则:

$ \dot V(x) \le - \left( {1 - {a^2}} \right){\mathit{\boldsymbol{X}}^{\rm{T}}}\mathit{\boldsymbol{X}} + \frac{{{\rho ^2}}}{2} $

(16)

从式(16)可以看出，只要系统的状态变量‖X‖_∞>ρ/$\sqrt {2\left( {1 - {a^2}} \right)} $，则$\dot V(x) < 0$成立，推导出整个航向控制系统(10)的一致最终有界。

为进一步讨论本文所提出的控制器对外界干扰的鲁棒性能，定义L₂增益鲁棒性能指标进行评价:

$ \int_0^t {{{\left\| \mathit{\boldsymbol{X}} \right\|}^2}} {\rm{d}}t \le {\mu _1}\int_0^t {{\mathit{\boldsymbol{w}}^2}(t){\rm{d}}t} + {\mu _2} $

(17)

式中μ₁和μ₂为无穷小的正数，则得到如下定理。

定理:在假设1~3的条件下，采用闭环增益成形算法，针对如式(10)所示的船舶航向控制系统进行控制器设计，在控制器参数T₁满足0＜T₁＜1时，系统的状态变量‖X‖_∞>ρ/$\sqrt {2\left( {1 - {a^2}} \right)} $，可保证整个控制器的一致最终有界，并可获得和控制器参数T₁相关的L₂增益鲁棒性能指标。

证明:

对式(16)从t=0到t=t₀进行积分，可得到:

$ \begin{array}{l} V\left( {{t_0}} \right) + \left( {1 - {a^2}} \right)\int_0^{{t_0}} {{\mathit{\boldsymbol{X}}^{\rm{T}}}\mathit{\boldsymbol{X}}{\rm{d}}t} - \\ \;\;\;\;\;\;\;\;\;\frac{1}{2}\int_0^{{t_0}} {\mathit{\boldsymbol{w}}{{(t)}^{\rm{T}}}\mathit{\boldsymbol{w}}(t){\rm{d}}t} \le V(0) \end{array} $

(18)

根据V(t₀)的定义，式(18)则可推导得到:

$ \int_0^{{t_0}} {{{\left\| \mathit{\boldsymbol{X}} \right\|}^2}} {\rm{d}}t \le \frac{1}{{2\left( {1 - {a^2}} \right)}}\int_0^{{t_0}} {{\mathit{\boldsymbol{w}}^2}(t){\rm{d}}t} + \frac{1}{{1 - {a^2}}}\varepsilon $

(19)

对式(19)中a的取值进行讨论，可得到:

$ a = \left\{ {\begin{array}{*{20}{l}} {{T_1},}&{{T_0} \ge 1}\\ {\left( {{T_1}\left( {1 + {T_0}{T_1}} \right)} \right)/\left( {{T_1} + {T_0}} \right),}&{0 < {T_0} < 1} \end{array}} \right. $

(20)

为了保证式(16)的负定，需保证1-a²>0，则根据式(20)可推导出控制器参数的取值范围为0＜T₁＜1。

因此，根据L₂增益鲁棒性能指标的定义(17)以及式(19)的结果可知，在航向控制系统的状态变量满足‖X‖_∞>ρ/$\sqrt {2\left( {1 - {a^2}} \right)} $条件下，只要在0＜T₁＜1范围内适当的调节闭环增益成形控制律(6)中的控制器参数T₁，可使船舶航向控制系统(10)获得如式(17)所示的L₂增益鲁棒性能指标，并且L₂增益鲁棒性能指标为1/(2(1-a²))，与控制器调节参数T₁相关。

证毕。

4 仿真试验与结果分析

为了验证本文所提出的航向控制器的稳定性和鲁棒性，本文选取大连海事大学的“育鲲”轮作为本文控制算法研究的被控对象，其基本参数可参考文献[7]。仿真实验采用4阶龙格-库塔方法进行数字仿真，仿真条件尽可能与航海实践相符。同时实验中为了精确的评价仿真实验的效果，本文参考文献[13]定义了控制器的控制效果评价指标CPEI(control performance evaluation index):

$ {\rm{CPEI}} = \sum\limits_{i = 0}^n {\left[ {{\alpha _1}{{\left( {\Delta {\psi _i}} \right)}^2} + {\alpha _2}\delta _i^2} \right]} $

式中:Δψ_i表示命令航向与当前航向的差值; δ_i表示舵角; 令α₁=10，α₂=0.1。

作为航向控制系统执行机构的舵机伺服系统是一个具有纯延迟、死区、滞环、饱和等非线性特性的电动液压系统，其对船舶航向控制的效果具有重要影响。因此为了保证本文仿真结果更接近航海实际情况，把舵机伺服系统近似成一阶惯性环节形式δ/δ_r=1/(T_rs+1)，并主要考虑其舵角限制和舵速限制，最终的仿真构建如图 2所示。

图 2中的舵机伺服系统主要由舵角限制器和舵速限制器2个部分组成，针对“育鲲”轮设为δ_max=30°，$\dot \delta $_max=2.8(°)/s。

假设船舶无外界干扰条件下进行仿真研究，初始航向设定为0°，目标航向为40°，仿真试验时间设定为800 s。则在此仿真场景下，为对本文所提出的控制器参数T₁选取一个合理的值，使船舶达到较好的航向控制效果。将控制器参数T₁在0＜T₁＜1的范围内按0.01的步长进行仿真实验，并计算每组CPEI值，得到如图 3所示的实验结果，进行数值回归取得CPEI为最小值时，控制器参数T₁为0.06，即控制器可在此处可获得较好的航向控制效果。

4.1 仿真1:大型船舶仿真实验

在海上航行时，高频海浪对船舶航向控制系统的输入干扰是引起船舶偏航的主要原因。其中三维不规则波可从海浪本质和机理的角度进行船舶航向控制系统进行仿真模拟，因此采用如下所示的合田改进JONSWAP波浪频谱:

$ S(f) = {\beta _j}h_{1/3}^2T_p^{ - 4}{f^{ - 5}}\exp \left[ { - \frac{5}{4}{{\left( {{T_p}f} \right)}^{ - 4}}} \right]{\gamma ^{\exp \left[ { - {{\left( {f/{f_p} - 1} \right)}^2}/2{\sigma ^2}} \right]}} $

(21)

式中γ为频谱参数，计算时取γ=3.3，

$ \begin{array}{l} {\beta _j} = \frac{{0.06238}}{{0.23 + 0.0336\gamma - 0.185{{(1.9 + \gamma )}^{ - 1}}}} \cdot \\ \;\;\;\;\;\;\;\;\left( {1.094 - 0.01915\ln \gamma } \right) \end{array} $

$ {T_p} = \frac{{{T_w}}}{{1 - 0.132{{(\gamma + 0.2)}^{ - 0.559}}}} $

$ \sigma = \left\{ {\begin{array}{*{20}{l}} {0.07}&{{f_w} \le {f_p}}\\ {0.09}&{{f_w} > {f_p}} \end{array}} \right. $

海浪方向函数如下:

$ G\left( {f,\theta } \right) = {\left[ {\int_{ - \frac{\pi }{2}}^{\frac{\pi }{2}} {{{\cos }^{2s}}\left( {\frac{\theta }{2}} \right){\rm{d}}\theta } } \right]^{ - 1}} \cdot {\left| {\cos \left( {\frac{{\theta - {\theta _0}}}{2}} \right)} \right|^{2s}} $

(22)

采用双叠加法模拟多向不规则波，即把多向不规则波的波能同时分布在一定频域和方向范围内，并把频域范围分成M份，方向分割成I份，同时把每个单元的波看成简谐波，其振幅表达式为:

$ {a_{mi}} = \sqrt {2s\left( {{\omega _m},{\theta _i}} \right)\Delta \omega \Delta \theta } $

其中，Δθ=(θ_max-θ_min)/I，Δω=(ω_H-ω_L)/M，最后进行多向不规则波的数值计算有:

$ \begin{array}{l} \eta \left( {x,y,t} \right) = \sum\limits_{m = 1}^M {\sum\limits_{i = 1}^I {{a_{mi}}} } \cos \left[ {{k_{mi}}\left( {x\cos {\theta _i} + } \right.} \right.\\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\left. {\left. {\gamma \sin {\theta _i}} \right) - {\omega _{mi}}t + {\varepsilon _{mi}}} \right] \end{array} $

式中ω_mi为波浪的角频率:

$ {\omega _{mi}} = {\omega _m} - \frac{1}{2}\Delta \omega + \left( {i - 1 + {\rm{RAN}}{{\rm{D}}_{mi}}} \right)\Delta \omega /I $

(23)

式中RAND_mi为(0, 1)范围内的均匀分布随机数。

假定仿真时平均海浪周期T_w=8 s，海浪有义波高h_1/3=3 m，则可得到如图 4所示的多向不规则波浪面仿真结果，主浪向为0°，坐标轴的原点为船舶重心。

本文设定船舶的初始航向为0°，目标航向为40°，根据式:

$ {I_{zz}}\dot r = \sum\limits_{m = 1,i = 1}^{M,I} {\frac{1}{{24}}{\rho _u}gBL\left( {{L^2} - {B^2}} \right)\sin (2\psi )a_{mi}^2(t)} $

(24)

和图 5所示的“育鲲轮”型线图计算船舶转艏向对波浪的响应，则得到船舶所受到的航向干扰，如图 6所示。式中:B、L分别表示为船舶的船宽、船长; ρ_w表示海水密度。

此时若控制器参数仍选为T₁=0.06，则得到如图 7所示的仿真结果。仿真结果与文献[14]非线性反馈控制方法进行比较。可看出本文在200 s后实现了良好的航向跟踪效果，且控制效果与文献[14]相差不大，CPEI值分别是5 785和6 498，本文在控制性能上约有12.32%的提升。

4.2 仿真2:小型无人船仿真实验

为验证本文所提出方法不仅适应于大型船舶，同样适用于小型无人船，因此利用文献[14]中的小型无人船“蓝信”号实验数据进行辨识所得到的模型数据进行仿真实验。其中“蓝信”号的运动数学模型如下:

$ {T_0}\dot r + r + a{r^3} = {K_0}\delta $

(25)

式中K₀=0.707，T₀=0.332，a=0.001。“蓝信”号的舵机模型采用如下:

$ \frac{\delta }{{{\delta _r}}} = \frac{{0.8955}}{{{s^2} + 1.5977s + 0.9702}} $

(26)

外界同样采用文献[14]中所示的五级风浪模型作为外界干扰，即白噪声驱动如下所示的二阶震荡传递函数来代表:

$ L(s) = \frac{{0.42s}}{{{s^2} + 0.36s + 0.37}} $

(27)

同时将式(25)中的非线性项ar³进行忽略，并利用本文的闭环增益成形算法进行控制器设计，其中控制器参数仍然取T₁=0.06，同时与文献[14]中的非线性修饰控制器进行对比仿真研究，仿真结果如图 8所示，2种方法所得到的CPEI值分别为17 538和23 509，本文控制方法在控制效果上会有13.04%的提升。

增大外界干扰到7级风浪，得到如图 9所示的控制仿真结果。同样控制器参数取T₁=0.06，并与文献[14]中的非线性修饰控制器进行对比仿真研究，2种方法所得到的CPEI值分别为63 509和66 165，本文控制方法在控制效果上有4.18%的效果提升，这主要是因为在恶劣海况下，小型船舶相比大型船舶对外界干扰会更敏感，鲁棒性会较大程度的变差。本文方法航向跟踪效果更好，但舵机消耗会增大。而文献[14]所采用的非线性修饰方法，在航向跟踪效果上会变差，但舵机消耗会变小。从控制器结构，计算复杂度等综合控制性能上分析，本文方法仍有一定优势。

5 结论

1) 在规划控制策略下采用闭环增益成形算法对无人船通用航向控制器进行设计，并采用李雅普诺夫稳定性理论和L₂增益鲁棒性能指标对控制器参数的范围进行了理论分析和证明。

2) 通过和非线性反馈控制方法进行仿真对比实验证明了本文所提出的航向保持控制器适用于不同尺度的无人船，特别是针对大型船舶本文方法约有12%的性能提升，而针对小型船舶本文方法越有4%的性能提升，而且本文控制方法结构简捷高效，在未来无人船舶规划控制策略中具有非常好的应用潜力。

未来工作要将本文所提出的航向反馈控制算法在实际无人船规划控制策略中加以应用，证明本文所提出方法的有效性和普适性。