0 引　言

非自航设备（例如受损的潜艇）的海上拖带运输方式主要有直接拖带、半潜船拖带和浮船坞拖带3种^[1]，如图1所示。直接拖带对非自航设备水面状态的密性和稳性要求较高，拖带过程存在一定的安全风险，而半潜船拖带吃水较大，对港区航道条件要求较高。浮船坞拖带可兼顾一定的稳性和不同的航道条件，适用面更广。海上拖带过程中，浮船坞在波浪作用下产生运动响应，但非自航设备特殊的圆形截面易导致其在运动过程中产生倾倒和滑移，影响安全，因此针对浮船坞拖带过程进行波浪响应分析，并设置可靠的系固结构是浮船坞拖带非自航设备时应考虑的关键因素。

现有针对浮船坞波浪响应的研究主要集中于系泊状态，董明海等^{[2 - 3]}、冯波^[4]等基于势流理论和悬链线理论就系泊浮船坞在极端海况下的运动响应和系泊揽拉力进行了分析。而针对浮船坞拖航的研究主要关注于空坞状态，李泽华等^[5]、胡维^[6]基于拖航实例，对空坞拖航作业过程中的重难点进行分析。同时拖航货物的系固结构研究主要围绕半潜船开展，杨宗俊等^[7]、刘子聃等^[8]对半潜船拖带重大件货物的系固结构抗滑移和翻转能力进行分析。但针对浮船坞装载着圆形截面的特殊货物进行海上拖带的研究还相对较少，有必要对拖带过程中的运动响应荷载及系固结构进行进一步研究。

本文采用非线性时域势流理论研究了装载SUBOFF潜艇的浮船坞在海上的运动响应，基于极端时刻的运动状态设计了货物系固结构，并基于水动力-结构有限元方法，将极端时刻的波浪及加速度荷载传递给结构有限元模型，分析了横向荷载对系固结构及浮船坞应力状态的影响，并针对性地进行了结构加强。本文研究可为圆形截面非自航设备的海上坞拖作业提供一定的参考。

1 计算原理与方法 1.1 运动响应短期预报

运动响应短期预报主要基于三维线性频域势流理论^[9]，假定海水为均匀、无粘不可压缩的理想流体。根据刚体动力学方程，运输船的运动方程可以写为：

$ \left[{\boldsymbol{M}}+{{\boldsymbol{M}}}_{a}(\omega )\right]\ddot{\eta }+{\boldsymbol{C}}\left(\omega \right)\dot{\eta }+{\boldsymbol{K}}\eta ={\boldsymbol{F}}\left(\omega \right)+{{\boldsymbol{F}}}_{T}\left(\omega \right)。$

(1)

式中：F(ω)为波浪引起的波激力；F_T(ω)拖船的拉力；M为船舶质量矩阵；M_a为附加质量矩阵；C为阻尼矩阵；K为刚度矩阵；$ \eta $为六自由度位移；ω为波浪频率。

假定波浪及船舶在波浪上的响应是线性的，引入速度势函数来描述流体的空间速度，并将速度势Φ分解为入射势Φ_I和扰动势Φ_P，系统满足理想流体无黏无旋不可压缩的拉普拉斯方程，以及线性自由面条件、物面条件、海底不可穿透条件。使用频域脉动格林函数，利用Hess-Smith均值面元法求解源强，通过边界积分即可获得速度势。得到速度势后，应用线性伯努利方程获得压力，并沿水动力表面进行压力积分，即可获得流体对运输船的1阶水动力，求解运动方程后即可得到幅值响应算子RAO：

$ {RAO}_{j}\left(\omega \right)=\frac{{x}_{j}(\omega )}{{\xi }_{a}(\omega )}。$

(2)

采用海况条件为Jonswap谱时，其谱密度表达式为：

$ {S}_{\xi }\left(\omega \right)=155\frac{H_{1/3}^{2}}{{{{T}_{1}}}^{4}{\omega }^{5}}{{\mathrm{e}}}^{\frac{-944}{{{{T}_{1}}}^{4}{\omega }^{4}}}{\gamma }^{{{e}^{-{{\left(\frac{0.191\omega {T}_{1}-1}{\sqrt{2}\sigma }\right)}^{2}}}}} 。$

(3)

式中：H_1/3为有义波高；T₁为波浪特征周期；$ \gamma $为谱峰因子；$ \sigma $=0.07（ω≤5.24/ T₁）或$ \sigma $=0.09（$ \omega $＞5.24/ T₁）。根据运动响应RAO和波浪谱按下式计算其响应谱$ {S}_{j}\left(\omega \right) $：

$ {S}_{j}\left(\omega \right)=RAO_{j}^{2}\left(\omega \right){S}_{\xi }\left(\omega \right) 。$

(4)

而响应谱的n阶谱距m_n：

$ {m}_{n}=\int _{0}^{\infty }{\omega }^{n}{S}_{j}\left(\omega \right)\mathrm{d}\omega。$

(5)

通过求得m₀、m₂及m₄，并进行谱宽修正后即可得到短期预报的有义值和极值。

1.2 非线性时域势流理论

三维线性频域势流理论无法考虑实际波面起浮，针对吃水浅、干舷低的浮船坞存在一定的偏差，而通过转变为时域计算使运动方程加入非线性Froud-Krylov项变为可能^[10]。运动方程引入脉冲响应函数（IRF）后可以写为：

$ \begin{split}&\left[{\boldsymbol{M}} + {\boldsymbol{M}}_a(\infty)\right]\ddot{{\boldsymbol{\eta}}}(t) + {\boldsymbol{C}}\dot{{\boldsymbol{\eta}}}(t) + \int_{0}^{t} {\boldsymbol{h}}(t-\tau)\dot{{\boldsymbol{\eta}}}(\tau){\rm{d}}\tau + \\&\qquad {\boldsymbol{K}}{{\eta}}(t) = {\boldsymbol{F}}(t) + {\boldsymbol{Q}}(t) + {\boldsymbol{F}}_t(t)。\end{split} $

(6)

式中：M_a(∞)为无穷大频率的流体附加质量矩阵；t为时间；τ为时间积分量；h为加速度脉冲响应函数矩阵；F(t)为船体所受一阶波浪力；Q(t)为船体所受二阶波浪力；F_T(t)为拖船的线性拉力。

考虑波浪在运输船水线附近的起伏变化，运动空间点坐标$ \overrightarrow{{X}} $上的压力P($ \overrightarrow{X} $,t)可以写波幅的函数：

$ P\left(\overrightarrow{X},t\right) = \text{Re}\left(\frac{i\omega {\xi }_{a}}{g}\varphi \left(\overrightarrow{X}\right){{\mathrm{e}}}^{-i\omega t}\right) + \text{Re}\left(-{\xi }_{a}\eta \left(\overrightarrow{X}\right){{\mathrm{e}}}^{-i\omega t}\right)。$

(7)

由于波浪起伏变化主要在高度方向，因此将某点$ P(\overrightarrow{X},t) $表达为P^z，将运输船和波浪间相对波幅定义为$ \xi (t) $，并将运输船静水线高度设为z=0，因此静水线处节点压力P⁰相对于静水时存在增大和减小两者情况，当波峰到达截面时P⁰＞0，修正后总压力可以表达为：

$ P\left(z,t\right)=\begin{cases} \rho g\xi (t)-\rho gz，{z}> 0，\\ \rho g\xi (t){{\mathrm{e}}}^{kz}-\rho gz，z< 0。\\ \end{cases} $

(8)

当波波到达截面时P⁰＜0，修正后总压力可以表达为：

$ P\left(z,t\right)=\begin{cases} 0，-\xi <{z}<0，\\ \rho g\xi (t){{\mathrm{e}}}^{kz}-\rho gz，{z}<-\xi 。\\ \end{cases} $

(9)

1.3 水动力-结构有限元方法

浮船坞的水动力荷载包括惯性力和水动压力，惯性力可在有限元计算中直接加载，而水动压力通过在有限元模型中船体湿表面网格点上进行直接计算，即在有限元模型中将湿表面网格进行“分离”，并在瞬态波浪面下湿表面网格上进行流体势函数计算，并通过伯努利方程进行积分，即可将水动力荷载“传递”给结构模型，进而进行结构有限元分析。

1.4 分析过程

本文基于BV船级社综合水动力-结构分析软件Hydrostar和HOMER，具体分析过程如图2所示。

2 运动响应预报 2.1 模型参数

浮船坞总长110.4 m，宽29 m，设计举力为50000 kN，浮船坞拖带2500 t SUBOFF潜艇缩比模型（以下简称货物）时吃水为2.68 m，排水量为8806 t，坞拖过程中详细参数见图3。

2.2 运动响应RAO

频域运动响应计算所采用的水动力网格如图4所示。

计算中的浪向定于如图5所示。

由于浮船坞拖航时航速较低，同时货物纵向采用钢缆固定，纵向系固难度较小，因此计算中将不考虑航速的影响，并主要关注于横荡、垂荡及横摇运动。浮船坞重心处运动及加速度响应RAO（幅频响应算子）计算结果如图6所示，其中横摇阻尼设置为5%。

可知，浮船坞在90°浪向下波浪频率为1.0 rad/s时，横摇幅值、垂荡及横摇加速度RAO曲线达到最大，该频率与浮船坞横摇固有频率0.97 rad/s对应，由此说明这一浪向和波浪频率对浮船坞横摇和垂荡运动的激励作用明显，是一个需要关注的危险工况。而浮船坞在90°浪向下横荡加速度RAO曲线呈现双峰，第一个峰值出现在0.84 rad/s处，此时波浪长度较长，对浮船坞横荡激励作用明显，因此也是一个需要关注的危险工况，而第二个峰值出现在1.54 rad/s处，此时波浪频率较高但波长较短，对浮船坞主要为持续的漫漂力作用，但对货物系固的影响较小，因此对该工况不进行进一步分析。

2.3 短期预报

将上述浮船坞运动响应RAO，基于JONSWAP波浪谱进行运动响应的短期预报。根据运输沿线海域资料^[11]，有义波高H_1/3为4.96 m，波浪特征周期T₁为7.81～18.2 s，γ值为1.0～5.0。将RAO结果带入波浪谱计算公式，并将3 h内概率为1/10的预报极值结果汇总形成表1。

3 时域水动力分析

通过频域计算可知，拖航过程中存在2个需关注的危险波浪频率，同时由于横荡加速度和横摇幅值对货物横向系固结构的影响最大，因此将针对上述危险波浪频率下，浮船坞的横荡加速度和横摇幅值进行进一步时域分析。

3.1 计算工况设置

应用等效设计波（RDW）^[12]计算法，即通过下式获得设计波的波幅A_RDW：

$ {A}_{{\rm{RDW}}}=\frac{X}{{\left| RAO({\omega }_{{\rm{RDW}}}，{\beta }_{{\rm{RDW}}})\right| }_{\max }}。$

(10)

式中：X为短期预报极值；$ {\left| RAO({\omega }_{{\rm{RDW}}}，{\beta }_{{\rm{RDW}}})\right| }_{\max } $为幅频响应算子RAO的最大值；$ {\omega }_{{\rm{RDW}}} $为等效设计波频率；$ {\beta }_{{\rm{RDW}}} $为等效设计波浪向。通过将横摇幅值、横荡加速度短期预报结果代入上式，得到如表2所示的2种工况设计波。

3.2 计算模型

浮船坞湿表面网格如图7所示。由于货物与浮船坞之间相对运动较整体的刚体运动为小量，因此计算时整体重心保持不变 。

3.3 计算结果

浮船坞重心处时域波浪运动响应曲线如图8所示，浮船坞在极端时刻波浪运动如图9所示。

图8(a)横摇角曲线中，工况1和工况2曲线的最大值基本相同，并与短期预报范围重合，而图8(b)横摇加速度曲线中，工况1曲线最大值18.8 °/s²明显大于工况2曲线最大值14.3 °/s，考虑到工况1中波浪的波幅和周期均小于工况2，由此说明工况1波浪频率对横摇运动的激励作用较工况2明显且与预报相符。

而图8(c)横荡加速度曲线可发现，工况1和工况2曲线的最大值均大于短期预报值1.42 m/s²，同时工况1曲线的最大值3.59 m/s²还大于工况2曲线的最大值2.35 m/s²。考虑到工况2频域计算中横荡加速度RAO大于工况1，且工况2波幅也大于工况1，由此说明采用非线性势流理论后，实际波面起伏对浮船坞横荡加速度的影响较大，基于线性频域计算的短期预报存在预报远小于实际横荡加速度的风险，同时还可说明工况1波浪频率对横荡加速度的激励作用强于工况2，是一个极为不利的危险工况。

综合图8中工况1的4个曲线，在横摇角达到最大19°时（5.31 s），横荡和横摇加速度也达到最大，但垂荡加速度较小。而垂荡加速度最大时（6.61 s），横摇角仅为1.7°，横荡和横摇加速度均较小。由此说明当浮船坞横摇角度达到最大时，横荡和横摇加速度也一同达到最大，此时货物所受的横向惯性力存在的协同叠加效应。由于工况1波浪频率与浮船坞横摇固有频域一致，综合上述分析可以说明，当拖航过程中遭遇与浮船坞固有频率相同的波浪且运动至横摇角度最大时刻，货物所受的横向耦合惯性力达到最大，是最可能产生倾倒和滑移的时刻，应针对该时刻的惯性力进行进一步的横向系固结构设计和校核。

此外，图8(d)垂荡加速度曲线中，工况1在正向半周期内曲线呈现双峰，而工况2曲线在正负半周期内均为单峰，结合图9横摇最大角度时刻的波浪运动图，工况1在最大横摇角时船体一侧底部出现露出水面、另一侧出现甲板上浪的现象，可以发现浮船坞在工况1时出现了底部拍击和甲板上浪现象，考虑到浮船坞本身具有吃水浅、前后干舷低的特点，由此说明虽然浮船坞拖带运输方式具有浅吃水适用面广的优点，但在特定波浪频率范围内容易出现底部拍击和甲板上浪现象，对浮船坞船体结构强度具有一定的影响。

4 结构分析 4.1 横向系固结构

由3.3节分析结果，当浮船坞横浪遭遇与其横摇固有频率相同的波浪，横摇角运动至最大时刻是货物最易产生倾倒和滑移的不利时刻，此时货物重心处横向加速度如表3所示。

通常浮船坞中坞墩结构按规范^[13]设置，如图10所示，由中墩和边墩组成。中墩起主支撑作用，边墩起辅支撑和侧扶作用，同时规范中仅对中墩承载能力进行了要求，对边墩承载能力未进行明确。由于拖航浮船坞最大横摇角度达到19°，且最大横向加速度达到了4.41 m/s²，常规边墩无法满足货物抗倾倒和滑移要求。因此本文在常规边墩基础上设置了横向系固结构，该结构由H型钢立柱、斜撑及肘板等结构组成，支撑点超过货物重心，可有效抵抗货物倾倒和滑移趋势。

4.2 结构计算模型

为验证横向系固结构的有效性，将3.3节中工况1最大横摇角时刻惯性荷载和波浪荷载传给有限元模型，进行结构有限元分析。浮船坞材料为B级钢，有限元模型采用壳单元和梁单元建立，网格大小为50 mm×50 mm。将货物在重心处简化为质量节点，同时考虑最大横摇时货物仅与一侧边墩及系固结构接触，将质量节点采用MPC约束与接触结构连接。结构有限元模型（隐去左舷甲板及舷墙）如图11所示。

4.3 计算结果

浮船坞合成应力云图如图12所示，计算结果汇总表如表4所示，其中K为材料许用系数^[14]。

通过图12(a)整体合成应力云图可以发现，浮船坞出现明显的扭转，但除支墩区域外，舷墙、外板及强框架处整体应力较小，说明结构整体满足承载能力需求。根据图12(b)局部应力云图可以发现，横向系固结构中应力分布合理，最大合成应力较小，说明系固加强结构满足潜艇横向抗倾倒及滑移需求。但一侧中墩和边墩结构甲板下横梁及纵桁存在较大应力区域，横梁腹板及纵行最大合成应力达到了268 MPa，横梁面板最大合成应力达到了192 MPa，超过了材料的许用应力，这是由于边墩不仅承受货物横向荷载，还承受垂向荷载，其所受总荷载值大于横向系固结构。

为保证拖带安全，有必要对浮船坞边墩下横梁、纵桁进行局部加强，将横梁腹板和纵桁加厚，横梁面板规格增大，并在边墩下纵桁处增设垂向T型加强筋，加强后的计算结果如图13和表4所示。通过云图及校核表可以说明，加强后浮船坞的局部高应力明显减小，其最大应力值满足材料许用要求且存在一定的余量，可以确保拖带时的结构安全。

5 结　语

本文针对装载非自航设备的浮船坞进行海上拖带的过程，进行了运动响应频域计算和短期预报，并对典型工况进行了非线性时域运动响应计算，并将典型时刻的惯性荷载和波浪荷载传递给有限元模型，计算了船坞、坞墩及系固结构的应力状态，通过分析得到如下结论：

1）当坞拖过程中横浪遭遇与浮船坞横摇固有频率接近的波浪时，浮船坞的横摇、横荡运动响应达到最大，最大横摇角、最大横荡加速度、最大横摇加速度在同一时刻出现，具有的耦合叠加效，此时波浪对货物的倾倒和滑移激励作用达到最大，应针对该时刻进行系固结构的设计和校核。

2）由于浮船坞具有浅吃水、前后低干舷的特点，横浪中大角度摇荡易造成底部拍击和甲板上浪现象，影响浮船坞的结构安全，实际作业过程中若长期横浪遭遇与浮船坞横摇固有频率接近的波浪环境时，应及时采用改变航向、调整吃水等措施进行规避。

3）通过设置本文所示的横向系固结构可减小货物出现倾倒和横向滑移的风险，但中墩及边墩仍然是主要支撑结构，作业前应根据浮船坞结构特点和其在目标海域的运动响应计算结果对甲板横梁、纵桁等局部部位进行加强，确保拖航安全。