0 引　言

随着海洋资源开发和深海科学研究的不断深入，水下机器人的应用日益广泛。自主遥控水下机器人（Autonomous and Remotely-operated Vehicle，ARV）结合了自主水下机器人（Autonomous Underwater Vehicle，AUV）和遥控水下机器人（Remotely operated Vehicle，ROV）的优点，具备在复杂海洋环境中灵活作业的能力。ARV被广泛应用于深海科学探测、海洋资源开发和环境监测，其利用微细光缆实现与母船之间的高速数据传输和实时控制，是完成复杂水下作业的重要保障^[1]。然而，目前ARV及其他使用微细光缆的水下装备主要采用抛弃式的微细光缆管理方式，在其完成任务后，将大量的微细光缆抛弃在海水中，可能会对海洋环境以及其他水下装备的作业产生一些未知的影响，所以开展微细光缆的回收及其相关技术问题的研究迫在眉睫。

目前，针对水下缆绳的动力学建模问题，主要方法包括有限差分法、有限段法和集中质量法等。Ablow等^[2]采用有限差分法，利用空间变量和时间变量组成动态微分方程，并采用隐式有限差分格式的数值解法。该模型允许使用较大的时间步长，运算速度较快。然而，针对节点数目及分段长度变化问题，有限差分法的灵活性较差。Xu等^{[3 − 4]}采用有限段法，建立了水下柔性段模型，将缆绳离散为一系列柔性弯曲段。通过对柔性段进行挠曲特性和外力分析，推导出力学方程。该模型考虑了缆绳的弯曲特性，精度较高，但由于计算量较大，难以应用于长距离缆绳的动力学仿真；Wang^[5]通过引入辅助节点，采用集中质量法，建立了拖缆回收和布放过程的动力学模型。该模型对恒定速度回收和布放缆绳的仿真效果较好，但由于分段水动力是通过相邻节点的平均速度计算的，导致在恒张力回收的启动阶段，仿真的数值波动较大。然而集中质量法具有良好的实时性，能够有效反映水下缆绳运动的瞬态特征。同时，该方法对于节点数目变化的适应性较强，处理起来较为灵活。

本文针对微细光缆动态回收过程，提出采用集中质量法建立微细光缆离散化模型，基于离散化模型对微细光缆进行了力学分析，建立微细光缆的动力学方程，并推演了微细光缆回收初始状态的力学方程以及计算方法。构建了基于微细光缆张力变化的节点控制策略，采用四阶龙格-库塔法对微细光缆的回收过程进行动力学仿真，揭示了微细光缆在回收过程中速度以及张力分布的动态变化规律，为水下机器人微细光缆管理优化和回收装置的设计提供了理论依据和前提条件。

1 微细光缆水下力学模型

在复杂多变的海洋环境下，长距离、大深度的水下微细光缆容易出现大幅度的弹性形变。由于其变形和受力过程涉及复杂的非线性特性，精确求解难度较大，因此需要通过简化模型，并结合逐步迭代的方法进行计算分析，以获取较为准确的结果^[6]。为此提出以下假设条件：

1）忽略由微细光缆内外层材质不同，引起的变形不均匀。

2）忽略微细光缆的弯曲刚度以及微细光缆的扭转。

3）假设微细光缆分段上的速度线性分布。

4）水动力对微细光缆的作用符合独立性原理，忽略切向水动力和法向水动力的耦合作用。

微细光缆的直径通常小于1 mm，弯曲刚度通常低于1×10⁻⁴ N·m²。根据弯矩-曲率方程，光缆弯曲产生的力矩等于光缆的弯曲刚度与曲率的乘积。这意味着在较大的回收张力和水动力作用下，微细光缆因弯曲产生的内力非常小。因此，为了简化模型计算，本文忽略了弯曲刚度在回收过程中对微细光缆运动产生的影响。

1.1 坐标系建立

如图1所示，将微细光缆分成$ n $段，分别记为$ {S_1} $,$ {S_2} $,…,$ {S_n} $。将节点分别记为$ {N_1} $,$ {N_2} $,…,$ {N_{n + 1}} $。其中回收端节点$ {N_{n + 1}} $与压坠器上搭载的回收装置相连，末端节点$ {N_1} $在初始状态下与ARV相连。回收阶段开始时，ARV本体与微细光缆分离。全局坐标系$ X - Y - Z $设置在回收装置上，$ X $轴和$ Y $轴沿着水平方向，$ Z $轴竖直向下。每个分段$ {S_i} $都有一个局部坐标系$ {b_i} - {n_i} - {t_i} $。假设节点$ {N_i} $的坐标为$\left( {{x_i},\;{y_i},\;{z_i}} \right)$，形变后分段$ {S_i} $长度记为${L_i}$：

$ {L_i} = \sqrt {{{\left( {{x_i} - {x_{i + 1}}} \right)}^2} + {{\left( {{y_i} - {y_{i + 1}}} \right)}^2} + {{\left( {{z_i} - {z_{i + 1}}} \right)}^2}}。$

(1)

则分段$ {S_i} $的方向向量$ {{\boldsymbol{t}}_i} $可以表示为：

$ {{\boldsymbol{t}}_i} = {\left[ {{{\boldsymbol{t}}_{xi}}}\ \ {{{\boldsymbol{t}}_{yi}}}\ \ {{{\boldsymbol{t}}_{zi}}} \right]^{\rm{T}}} = \frac{1}{{{L_i}}}{\left[ {{x_i} - {x_{i + 1}}}\ \ {{y_i} - {y_{i + 1}}}\ \ {{z_i} - {z_{i + 1}}} \right]^{\rm{T}}}。$

(2)

由于忽略扭转，由全局坐标系$ X - Y - Z $到局部坐标系$ {b_i} - {n_i} - {t_i} $的旋转矩阵${R_i}$可以表示为：

$ {{\boldsymbol{R}}_i} = \left[ {\begin{array}{*{20}{c}} {\displaystyle\frac{{{t_{zi}}}}{{\sqrt {1 - t_{yi}^2} }}}&{\displaystyle\frac{{ - {t_{xi}}{t_{yi}}}}{{\sqrt {1 - t_{yi}^2} }}}&{{t_{xi}}} \\ 0&{\sqrt {1 - t_{yi}^2} }&{{t_{yi}}} \\ {\displaystyle\frac{{ - {t_{xi}}}}{{\sqrt {1 - t_{yi}^2} }}}&{\displaystyle\frac{{ - {t_{zi}}{t_{yi}}}}{{\sqrt {1 - t_{yi}^2} }}}&{{t_{zi}}} \end{array}} \right]，$

(3)

$ \left\{ {\begin{aligned} &{{b_i} = {{\left[ {\begin{array}{*{20}{c}} {\displaystyle\frac{{{t_{zi}}}}{{\sqrt {1 - t_{yi}^2} }}}\ \ \ 0\ \ \ {\displaystyle\frac{{ - {t_{xi}}}}{{\sqrt {1 - t_{yi}^2} }}} \end{array}} \right]}^{\rm{T}}}}，\\ &{{n_i} = {{\left[ {\begin{array}{*{20}{c}} {\displaystyle\frac{{ - {t_{xi}}{t_{yi}}}}{{\sqrt {1 - t_{yi}^2} }}}\ \ \ {\sqrt {1 - t_{yi}^2} }\ \ \ {\displaystyle\frac{{ - {t_{zi}}{t_{yi}}}}{{\sqrt {1 - t_{yi}^2} }}} \end{array}} \right]}^{\rm{T}}}} 。\end{aligned}} \right. $

(4)

海水中微细光缆在外力的作用下会产生弹性变形，微细光缆泊松比为0.5。未形变的微细光缆总长度记为$ S $，则分段$ {S_i} $未形变的长度${l_i}$为：

$ {l_i} = \frac{S}{n}。$

(5)

分段$ {S_i} $的应变，记为${\varepsilon _i}$：

$ {\varepsilon _i} = \frac{{{L_i} - {l_i}}}{{{l_i}}}。$

(6)

1.2 受力分析

微细光缆在海水中受到的外力主要包含回收端和水下机器人端的拉力、重力、浮力以及水动力^[7]，将这些外力集中到模型的质量点上，结合质量点之间的张力作用，对微细光缆进行力学分析。

$ E $为微细光缆的弹性模量，分段$ \mathit{{S}_{{i}}} $未形变的横截面积为$ {{a}} $，分段$ {S_i} $的张力${T_i}$为：

$ {T_i} = Ea{\varepsilon _i}。$

(7)

由于节点$ {N_i} $受到两侧分段的张力作用，则节点$ {N_i} $受到的张力作用力$ {{\boldsymbol{F}}_{T,i}} $为：

$ {{\boldsymbol{F}}_{T,i}} = {T_{i - 1}}{{\boldsymbol{t}}_{i - 1}} - {T_i}{{\boldsymbol{t}}_i}。$

(8)

将分段$ {{{S}}_{{i}}} $受到的重力和浮力平均分配到两端节点上，微细光缆单位长度的质量为$ m $，$ \rho $为海水的密度，$ g $为重力加速度，则节点$ {N_i} $受到的净重力$ \boldsymbol{F}_{G,i} $为：

$ {{\boldsymbol{F}}_{G,i}} = \frac{1}{2}g\left( {m - \rho a} \right)\left( {{l_i} + {l_{i - 1}}} \right){\left[ 0\ \ 0\ \ 1 \right]^{\rm{T}}}。$

(9)

微细光缆分段是细长圆柱体，可以使用半经验的莫里森公式计算水动力。微细光缆分段$ {S_i} $受到的水动力可以分为切向水动力和法向水动力。假设节点$ {N_i} $所在坐标的海水流动速度为${{\boldsymbol{V}}_{w,i}}$，节点$ {N_i} $的速度为${\dot {\boldsymbol{X}}_i}$。则在微细光缆分段$ {S_i} $三个方向上的计算水流速度为：

$ {\left[ {{V_{sb,i}}}\ \ \ {{V_{sn,i}}}\ \ \ {{V_{st,i}}} \right]^{\rm{T}}} = {\boldsymbol{R}}_i^{\rm{T}}({{\boldsymbol{V}}_{w,i}} - {\dot {\boldsymbol{X}}_i})。$

(10)

在计算分段受到的水动力时，通常使用分段两端节点的平均速度来计算整个分段的水动力，再将分段的水动力平均分配到两端节点上。这种方法可以应用于分段两端节点速度差异不大的场景^[8]。但对于采用恒张力回收微细光缆的场景，随着缆绳长度的减少，节点整体运动是加速运动，采用这种方法会使下端节点受到的水动力偏大，引起数值上的波动。本文通过引入速度比例因子，利用微积分的方法计算节点上分配的水动力。微细光缆分段$ {S_i} $的切向速度比例因子$ {k_{t,i}} $为：

$ {k}_{t,i}=\left\{\begin{aligned}&\displaystyle\frac{\left|{V}_{st,i+1}\right|}{\left|{V}_{st,i+1}\right|+\left|{V}_{st,i}\right|},{V}_{st,i+1}和{V}_{st,i}不同时为0，\\ &\displaystyle\frac{1}{2},{V}_{st,i+1}和{V}_{st,i}同时为0。\end{aligned}\right. $

(11)

$ {\rm{sign}}(x) = \left\{ \begin{aligned} &1 ,x \gt 0 ，\\ &0 ,x = 0 ，\\ &- 1 ,x \lt 0。\end{aligned} \right. $

(12)

在微细光缆分段$ {S_i} $的切线方向上，分段两端节点的水动力分别为：

$ {\left\{ \begin{aligned} & {{\boldsymbol{H}}_{t,i,i + 1}} = \displaystyle\frac{1}{2}{C_t}\rho {\text{π}} \displaystyle\frac{{{d_i}}}{{\sqrt {1 + {\varepsilon _i}} }}{\rm{sign}}({V_{st,i + 1}})\times\\ &\int_0^{{k_{t,i}}{L_i}} {\left( {V_{st,i + 1}^2 + 2\displaystyle\frac{{{V_{st,i + 1}}{V_{st,i}} - V_{st,i + 1}^2}}{{{L_i}}}x + \displaystyle\frac{{{{({V_{st,i}} - {V_{st,i + 1}})}^2}}}{{L_i^2}}{x^2}} \right)} {\rm{d}}x \cdot {{\boldsymbol{t}}_i} ，\\ &{{\boldsymbol{H}}_{t,i,i}} = \displaystyle\frac{1}{2}{C_t}\rho{\text{π}} \frac{{{d_i}}}{{\sqrt {1 + {\varepsilon _i}} }}{\rm{sign}}({V_{st,i}})\times\\ & \int_{{k_{t,i}}{L_i}}^{{L_i}} {\left( {V_{st,i + 1}^2 + 2\displaystyle\frac{{{V_{st,i + 1}}{V_{st,i}} - V_{st,i + 1}^2}}{{{L_i}}}x + \displaystyle\frac{{{{({V_{st,i}} - {V_{st,i + 1}})}^2}}}{{L_i^2}}{x^2}} \right)} {\rm{d}}x \cdot {{\boldsymbol{t}}_i}。\end{aligned} \right.} $

(13)

同理，在微细光缆分段$ {S_i} $的法线方向上，分段两端节点的水动力分别为：

$ k_{b,i}=\left\{\begin{aligned} & \displaystyle\frac{\left|V_{sb,i+1}\right|}{\left|V_{sb,i+1}\right|+\left|V_{sb,i}\right|},V_{sb,i+1}和V_{sb,i}不同时为0，\\ & \displaystyle\frac{1}{2},V_{sb,i+1}和V_{sb,i}同时为0。\end{aligned}\right. $

(14)

$ {\left\{ \begin{aligned} &{{\boldsymbol{H}}_{b,i,i + 1}} = \displaystyle\frac{1}{2}{C_n}\rho \displaystyle\frac{{{d_i}}}{{\sqrt {1 + {\varepsilon _i}} }}{\rm{sign}}({V_{sb,i + 1}})\times\\ &\int_0^{{k_{b,i}}{L_i}} {\left( {V_{sb,i + 1}^2 + 2\displaystyle\frac{{{V_{sb,i + 1}}{V_{sb,i}} - V_{sb,i + 1}^2}}{{{L_i}}}x + \displaystyle\frac{{{{({V_{sb,i}} - {V_{sb,i + 1}})}^2}}}{{L_i^2}}{x^2}} \right)} {\rm{d}}x \cdot {{\boldsymbol{b}}_i}，\\ &{{\boldsymbol{H}}_{b,i,i}} = \frac{1}{2}{C_n}\rho \displaystyle\frac{{{d_i}}}{{\sqrt {1 + {\varepsilon _i}} }}{\rm{sign}}({V_{sb,i}})\times\\ &\int_{{k_{b,i}}{L_i}}^{{L_i}} {\left( {V_{sb,i + 1}^2 + 2\displaystyle\frac{{{V_{sb,i + 1}}{V_{sb,i}} - V_{sb,i + 1}^2}}{{{L_i}}}x + \displaystyle\frac{{{{({V_{sb,i}} - {V_{sb,i + 1}})}^2}}}{{L_i^2}}{x^2}} \right)} {\rm{d}}x \cdot {{\boldsymbol{b}}_i} 。\end{aligned} \right.} $

(15)

$ {k}_{n,i}=\left\{\begin{aligned}&\displaystyle\frac{\left|{V}_{sn,i+1}\right|}{\left|{V}_{sn,i+1}\right|+\left|{V}_{sn,i}\right|},{V}_{sn,i+1}和{V}_{sn,i}不同时为0，\\ &\displaystyle\frac{1}{2},{V}_{sn,i+1}和{V}_{sn,i}同时为0。\end{aligned} \right.$

(16)

$ {\left\{ \begin{aligned} &{{\boldsymbol{H}}_{n,i,i + 1}} = \displaystyle\frac{1}{2}{C_n}\rho \displaystyle\frac{{{d_i}}}{{\sqrt {1 + {\varepsilon _i}} }}{\rm{sign}}({V_{sn,i + 1}})\times\\ &\int_0^{{k_{n,i}}{L_i}} {\left( {V_{sn,i + 1}^2 + 2\displaystyle\frac{{{V_{sn,i + 1}}{V_{sn,i}} - V_{sn,i + 1}^2}}{{{L_i}}}x + \displaystyle\frac{{{{({V_{sn,i}} - {V_{sn,i + 1}})}^2}}}{{L_i^2}}{x^2}} \right)} {\rm{d}}x \cdot {{\boldsymbol{n}}_i}，\\ &{{\boldsymbol{H}}_{n,i,i}} = \displaystyle\frac{1}{2}{C_n}\rho \displaystyle\frac{{{d_i}}}{{\sqrt {1 + {\varepsilon _i}} }}{\rm{sign}}({V_{sn,i}})\times\\ &\int_{{k_{b,i}}{L_i}}^{{L_i}} {\left( {V_{sn,i + 1}^2 + 2\displaystyle\frac{{{V_{sn,i + 1}}{V_{sn,i}} - V_{sn,i + 1}^2}}{{{L_i}}}x + \displaystyle\frac{{{{({V_{sn,i}} - {V_{sn,i + 1}})}^2}}}{{L_i^2}}{x^2}} \right)} {\rm{d}}x \cdot {{\boldsymbol{n}}_i} 。\end{aligned} \right. }$

(17)

节点$ {N_i} $受到的水动力$ {{\boldsymbol{F}}_{H,i}} $：

$ {{\boldsymbol{F}}_{H,i}} = {{\boldsymbol{H}}_{t,i,i}} + {{\boldsymbol{H}}_{n,i,i}} + {{\boldsymbol{H}}_{b,i,i}} + {{\boldsymbol{H}}_{t,i - 1,i}} + {{\boldsymbol{H}}_{n,i - 1,i}} + {{\boldsymbol{H}}_{b,i - 1,i}}。$

(18)

节点${N_i}$受到的合外力${{\boldsymbol{F}}_i}$：

$ {{\boldsymbol{F}}_i} = {{\boldsymbol{F}}_{T,i}} + {{\boldsymbol{F}}_{G,i}} + {{\boldsymbol{F}}_{H,i}}。$

(19)

对于回收端节点$ {N_{n + 1}} $与末端节点$ {N_1} $，只有一条分段与其连接，在建立力学方程时只需对与其相连的分段进行受力分析。

2 微细光缆回收的初始状态

微细光缆的回收过程是一个时变的动态过程。从微细光缆与ARV本体分离时刻起，即可启动其基于恒张力的回收，此时需要研究并分析微细光缆的初始状态。即动态回收开始时，微细光缆各节点坐标以及速度^[9]。本文是以微细光缆处于静止状态（各节点速度为0），作为微细光缆回收动力学仿真的初始状态。

2.1 初始状态的力学方程

微细光缆回收端节点$ {N_{n + 1}} $与回收装置连接，末端节点$ {N_1} $与ARV本体连接，这2个点固定，不需要对其展开受力分析，其余节点在净重力、张力、水动力的作用下，处于静止状态。将这些节点坐标写为矩阵的形式：

$ {\boldsymbol{X}} = {\left( {{x_2},{y_2},{z_2};{x_3},{y_3},{z_3}; \cdots ;{x_n},{y_n},{z_n}} \right)^{\rm{T}}}。$

(20)

微细光缆回收端节点$ {N_{n + 1}} $的坐标为$ （0,0,0） $，假设ARV本体的坐标（节点$ {N_1} $）为$ ({x_1},{y_1},{z_1}) $，所有节点的速度都为$ 0\;{\text{m/s}} $，由式(19)，求得力学方程为：

$ {\boldsymbol{F}}({\boldsymbol{X}}) = \left\{ {\begin{array}{*{20}{c}} {{{\boldsymbol{F}}_2}({\boldsymbol{X}}) = 0}，\\ {{{\boldsymbol{F}}_3}({\boldsymbol{X}}) = 0}，\\ \vdots \\ {{{\boldsymbol{F}}_n}({\boldsymbol{X}}) = 0}。\end{array}} \right. $

(21)

2.2 初始状态的求解

采用牛顿法求解微细光缆回收初始状态的力学方程。本文采用0.4 mm的微细光缆，微细光缆的物理参数如表1，假设水流大小为沿X轴−0.1 m/s，微细光缆长度为6000 m，ARV本体的起始坐标为(1000, 0, 4000)。

微细光缆的水下形态如图2(a)所示，在外力作用下，光缆呈现开口向右的弧形。由于水流速度较小，部分微细光缆的坐标深度超过4000 $ {\text{m}} $，这主要是重力作用的结果。张力分布曲线如图2(b)所示。微细光缆的两端分别固定在回收装置和ARV本体上，需要承受较大的张力以平衡作用在微细光缆上的外力。同时，微细光缆的大部分重力集中在回收端节点上。因此随着微细光缆距离回收点长度的增加，张力呈现出先减小后增大的趋势，而且在回收端和ARV端，张力相对较大。

3 微细光缆回收过程动力学建模

为了实现回收效率最大化，应选择接近微细光缆材料破断力的拉力值作为恒张力回收的张力值，确保微细光缆在回收过程中不会发生断裂^[10]。为此，通过仿真分析可以有效计算出回收速度，以此作为实际回收过程中的操作参考，提高回收过程的安全性与可靠性，同时最大限度地提升回收效率。

3.1 动力学模型构建与分析

物体在流体中运动，由物体加速度引起的流体作用力通常称为附加质量力。微细光缆分段是细长圆柱体，当其在海水中做加速运动时，切向上引起的附加质量力相比法向上非常小，可以忽略不计^[11]。假设分段$ {S_i} $加速度为${\ddot {\boldsymbol{L}}_{S,i}}$，附加质量系数为$ {k_s} $，附加质量力${{\boldsymbol{F}}_{S,i}}$可表示为：

$ {{\boldsymbol{F}}_{S,i}} = {{\boldsymbol{R}}_i}{k_s}\rho a{l_i}\left[ {\begin{array}{*{20}{c}} 1&0&0 \\ 0&1&0 \\ 0&0&0 \end{array}} \right]{{\boldsymbol{R}}_i}^{\rm{T}}{\ddot {\boldsymbol{L}}_{S,i}}。$

(22)

则分段$ {S_i} $的附加质量矩阵${{\boldsymbol{M}}_{as,i}}$为：

$ {{\boldsymbol{M}}_{as,i}} = \frac{1}{4}{\text{π}} {k_s}\rho {d^2}{l_i}{{\boldsymbol{R}}_i}\left[ {\begin{array}{*{20}{c}} 1&0&0 \\ 0&1&0 \\ 0&0&0 \end{array}} \right]{{\boldsymbol{R}}_i}^{\rm{T}}。$

(23)

节点${N_i}$的附加质量矩阵${{\boldsymbol{M}}_{a,i}}$可以表示为：

$ {{\boldsymbol{M}}_{a,i}} = \frac{1}{2}({{\boldsymbol{M}}_{as,i}} + {{\boldsymbol{M}}_{as,i - 1}})。$

(24)

则节点${N_i}$的质量矩阵${{\boldsymbol{M}}_{n,i}}$为：

$ {{\boldsymbol{M}}_{n,i}} = {{\boldsymbol{M}}_{a,i}} + \frac{1}{2}m({l_i} + {l_{i - 1}})\left[ {\begin{array}{*{20}{c}} 1&0&0 \\ 0&1&0 \\ 0&0&1 \end{array}} \right] 。$

(25)

节点${N_i}$在合外力${{\boldsymbol{F}}_i}$作用下产生的加速度${\ddot {\boldsymbol{X}}_i}$为：

$ {\ddot {\boldsymbol{X}}_i} = {\boldsymbol{M}}_{n,i}^ - {{\boldsymbol{F}}_i}。$

(26)

在微细光缆的回收过程中，对于已回收的微细光缆的受力不做考虑。随着海水中微细光缆长度的不断减少，对于微细光缆离散化模型，在不改变所有分段长度的前提下，必然引起节点数量的减少。所以在仿真的过程中需要动态的减少微细光缆的节点数量，同时动态调整以形成新的动力学状态方程。为了模拟微细光缆回收的动态过程，以上文中的仿真结果作为初始状态。回收开始时，ARV本体与微细光缆分离，同时在回收端节点上增加一个回收张力$ {F_0} $。在回收过程中，回收张力的方向随着分段$ {S_n} $切线方向动态变化，节点$ {N_{n + 1}} $受力方程为：

$\begin{split} {{\boldsymbol{M}}_{n,n + 1}}{\ddot {\boldsymbol{X}}_{n + 1}} =& {{\boldsymbol{F}}_{T,n + 1}} + {{\boldsymbol{F}}_{G,n + 1}} + {{\boldsymbol{F}}_{H,n + 1}} - {F_0}{{\boldsymbol{t}}_n} = \\&{T_n}{{\boldsymbol{t}}_n} - {F_0}{{\boldsymbol{t}}_n} + {{\boldsymbol{F}}_{G,n + 1}} + {{\boldsymbol{F}}_{H,n + 1}}。\end{split} $

(27)

节点$ {N_n} $受力方程为：

$ {{{\boldsymbol{M}}_{n,n}}{\ddot {\boldsymbol{X}}_n} = {{\boldsymbol{F}}_{T,n}} + {{\boldsymbol{F}}_{G,n}} + {{\boldsymbol{F}}_{H,n}} = {T_{n - 1}}{{\boldsymbol{t}}_{n - 1}} - {T_n}{{\boldsymbol{t}}_n} + {{\boldsymbol{F}}_{G,n}} + {{\boldsymbol{F}}_{H,n}}。} $

(28)

式中：${M_{n,n + 1}}$、${F_{G,n + 1}}$、${F_{H,n + 1}}$、${M_{n,n}}$均为长度函数，在分段$ {S_n} $的回收过程中，随分段$ {S_n} $的长度动态变化：

${ \left\{ \begin{array}{*{20}{l}} {{l_n} \to 0}，\\ {{M_{n,n + 1}}:\displaystyle\frac{1}{2}{M_{as,n}} + \displaystyle\frac{1}{2}m{l_n}\left[ {\begin{array}{*{20}{c}} 1&0&0 \\ 0&1&0 \\ 0&0&1 \end{array}} \right] \to 0} ，\\ {{F_{G,n + 1}} \to 0}，\\ {{F_{H,n + 1}} \to 0}，\\ {M_{n,n}}:\displaystyle\frac{1}{2}({M_{as,n}} + {M_{as,n - 1}}) + \displaystyle\frac{1}{2}m({l_n} + {l_{n - 1}})\left[ {\begin{array}{*{20}{c}} 1&0&0 \\ 0&1&0 \\ 0&0&1 \end{array}} \right] \to \\ \displaystyle\frac{1}{2}{M_{as,n - 1}} + \displaystyle\frac{1}{2}m{l_{n - 1}}\left[ {\begin{array}{*{20}{c}} 1&0&0 \\ 0&1&0 \\ 0&0&1 \end{array}} \right] 。\end{array} \right.} $

(29)

在分段$ {S_n} $的回收过程中，分段$ {S_n} $的张力$ {T_n} $是逐渐增大的，根据式(27)，结合上述物理量的变化，当分段$ {S_n} $长度趋近于0时，分段$ {S_n} $的张力$ {T_n} $趋近于回收张力$ {F_0} $。此时进行删除节点的操作，相当于重新在节点$ {N_n} $上施加一个沿着切线方向的回收张力，对微细光缆回收过程动力学仿真的影响较小。

3.2 动力学模型求解方法

基于微细光缆的初始状态，设置合理的时间步长，利用四阶龙格-库塔法，对一个时间步长之后微细光缆的状态进行求解。计算此时分段$ {S_n} $的长度$ {L_n} $：

$ {L_n} = \sqrt {{{\left( {{x_n} - {x_{n + 1}}} \right)}^2} + {{\left( {{y_n} - {y_{n + 1}}} \right)}^2} + {{\left( {{z_n} - {z_{n + 1}}} \right)}^2}}。$

(30)

微细光缆回收装置的坐标为$ (0,0,0) $。在模型中删除已经回收的微细光缆，将节点$ {N_{n + 1}} $坐标替换为$ (0,0,0) $，计算替换后分段$ {S_n} $的长度$ {L_{n,k}} $：

$ {L_{n,k}} = \sqrt {x_n^2 + y_n^2 + z_n^2}。$

(31)

假设替换前分段$ {S_n} $未形变的长度为$ {l_n} $。按照替换前后，分段$ {S_n} $张力不变原则，计算替换后分段$ {S_n} $未形变的长度$ {l_{n,k}} $以及分段$ {S_n} $的张力$ {T_{n,k}} $：

$ \left\{ {\begin{aligned} &{{l_{n,k}} = \displaystyle\frac{{{L_{n,k}}}}{{{L_n}}}{l_n}}，\\ &{{T_{n,k}} = \displaystyle\frac{{Ea({L_{n,k}} - {l_{n,k}})}}{{{l_{n,k}}}}} 。\end{aligned}} \right. $

(32)

替换前节点$ {N_{n + 1}} $的速度为${\dot {\boldsymbol{X}}_{n + 1}}$，节点$ {N_n} $的速度为${\dot {\boldsymbol{X}}_n}$。利用线性插值的方法，计算替换后节点$ {N_{n + 1}} $的速度${\dot {\boldsymbol{X}}_{n + 1,k}}$（即微细光缆的回收速度）：

$ {\dot {\boldsymbol{X}}_{n + 1,k}} = \frac{{{L_{n,k}}}}{{{L_n}}}({\dot {\boldsymbol{X}}_{n + 1}} - {\dot {\boldsymbol{X}}_n}) + {\dot {\boldsymbol{X}}_n}。$

(33)

通过理论模型分析，在微细光缆回收的过程中，分段$ {S_n} $的长度趋近于0时，分段$ {S_n} $的张力趋近于回收张力$ {F_0} $。但是实际上由于初始状态力学方程的解是数值解，当分段$ {S_n} $的长度减小到一定程度后，会放大数值解存在的误差影响，导致数值的不稳定。通常的方法是设置一个长度值$ {{{L}}_0} $，当分段$ {S_n} $的长度小于$ {{{L}}_0} $时，进行删除节点的操作。每次删除节点后，节点$ {N_{n + 1}} $受到张力作用力的跃升幅度不固定。本文提出基于分段$ {S_n} $张力变化的节点控制策略，通过需要设置张力阈值$ {T_0} $。当分段$ {S_n} $的张力>$ {T_0} $时，删除节点，保证每次删除节点后，节点$ {N_{n + 1}} $受到的张力作用力跃升幅度是一个较小且稳定的数值。将张力阈值$ {T_0} $设置为回收张力$ {F_0} $数值的99.95%以上。当分段$ {S_n} $的张力$ {T_{n,k}} $>$ {T_0} $时，删除现有的节点$ {N_n} $，计算删除后的节点数目$n$以及分段$ {S_n} $未形变的长度$ {l_{n,k}} $：

$ \left\{ {\begin{aligned} &{{l_{n,k}} = {l_{n,k}} + {l_{n - 1}}}，\\ &{n = n - 1} 。\end{aligned}} \right. $

(34)

重复上述过程，直到当节点数$n = 0$时，微细光缆回收过程动态仿真完成。

3.3 回收过程的动态仿真分析

基于上文所述微细光缆的初始状态，采用接近微细光缆材料破断力的拉力值（即100 $ {\text{N}} $）作为回收过程动力学仿真的张力值。仿真结果表明，微细光缆回收过程持续了2741.68 s。回收过程中微细光缆的形态变化如图3所示。在回收开始时，由于失去了ARV本体对末端的约束作用，末端节点受外力影响发生运动。随着光缆长度的减少，微细光缆的弧度逐渐减小，最终趋近于直线。由于海水流速较低，水动力对于微细光缆运动的影响小于回收张力的影响，同时大部分光缆位于回收点的右侧。所以在回收的前期，微细光缆的运动是逐渐向右侧偏移的。但在回收接近尾声时，微细光缆明显向左侧发生偏移，这就是水流作用的结果。

通过对微细光缆回收过程进行动力学仿真，得到微细光缆回收速度变化曲线如图4所示，在回收开始时，回收速度迅速上升，并在短时间内达到峰值。这是由于回收开始前，微细光缆内部的张力较小，在瞬时施加较大的回收张力后，微细光缆的弹性作用使得回收速度在短时间内迅速增加。随着回收的进行，回收速度逐渐减小，这是由于张力和水动力对回收端节点的运动产生了阻尼作用，从而导致回收速度减小。图4(c)为部分回收速度下降阶段的速度分布曲线，速度从回收端到末端递减。随着回收的进行，末端节点速度逐渐增加。靠近末端的速度曲线发生明显的弯曲，这是由于微细光缆在水流作用下呈弯曲状态，在回收张力作用下，微细光缆弯曲弧度逐渐减小，使得靠近末端节点速度变化相对较慢；随着时间的推移，回收速度停止下降并逐渐加速上升，最终理论速度趋近无穷大。这是由于微细光缆的长度减少到一定程度后，水动力作用下整体阻力减小的结果。图4(d)为部分回收速度上升阶段微细光缆的速度分布曲线，曲线间隔逐渐增加，表明微细光缆加速度的增加。在实际的回收过程中，采用接近光缆材料破断力的张力进行回收，可以保证回收效率的最大化，但是为了避免微细光缆发生断裂，实际回收速度不应超过仿真回收速度。

图5(a)为回收速度下降阶段的张力分布曲线，张力呈曲线分布。在施加回收张力后，运动逐渐向下传递，节点坐标的变化导致分段张力的变化，各节点所受张力作用力逐渐趋于一致。靠近末端，张力分布曲线的斜率变得平缓，这与末端节点速度变化较慢相对应。图5(b)为回收速度上升阶段的张力分布曲线，张力呈现线性分布，随着光缆长度的减少，斜率逐渐增大，说明节点所受张力作用力逐渐增大。通过回收速度变化曲线可以看出，这段时间内，回收速度曲线的斜率逐渐增加，说明回收运动的加速度在增大。节点速度的增加导致了分段水阻力的增大，因此只有节点受到张力作用力增大才能引起回收加速度的增大。仿真结果表明，微细光缆张力的最大值处于回收点位置，数值上等于施加的回收张力，所以在实际回收的过程中，回收张力应小于微细光缆的破断力。

4 结　语

1）微细光缆回收速度分布变化规律：回收开始时，由于光缆自身弹性作用，回收速度快速达到峰值后，逐渐下降至稳态值；随着微细光缆剩余长度的不断减少，回收速度开始逐渐上升并且加速度不断增大，最后理论的回收速度趋于无穷大；根据回收速度及其变化规律，可进一步推算出采用恒张力回收特定长度的微细光缆所需时间，从而为微细光缆回收效率的提升提供参考和指导。

2）微细光缆张力分布变化规律：回收的过程中，张力由非线性分布逐渐趋于线性分布，回收端节点受到的张力最大，保证回收端张力不超过光缆的破断力，即可保证微细光缆的安全回收，进而为水下机器人微细光缆管理能力的提升和回收装备的研发提供了理论参考和设计依据。