近年来，由于海洋渔业资源的减少以及人们对海产品需求的增加，海洋养殖业正快速发展。网箱养殖作为一种主要的养殖方式，在我国海水鱼类养殖中占有重要地位，2019年我国海水鱼类养殖产量为160.6×10⁴ t，其中网箱养殖产量占比47.1%^[1]。为了深入研究高海况下海洋荷载对网箱性能的影响，学者们基于数值模拟、物理模型试验和现场实测对网箱系统的水动力特性进行了大量分析，以往的研究涵盖了网箱的所有构件，包括浮架、网衣和系泊系统。Lee等^[2]通过数值模型为浮架结构设计提供了有价值的信息。Huang和Pan^[3]基于长周期环境载荷评估了单点系泊系统的锚绳失效风险，给出了使用期内的失效概率。Shainee等^[4]研究了一种自潜式单点系泊网箱系统在波流下的淹没特性。Fredriksson等^[5]利用数值模型和现场实测数据，研究了大型渔场在无波情况下的系泊系统张力。Zhao等^[6-7]基于集中质量法和刚体运动学原理对箱型网箱和圆柱形网箱在波浪和水流中的动力行为进行了分析。Decew等^[8]提出了一种声学方法，用于监测在水流中部署的小型网箱的运动和变形。

然而，上述大部分研究工作是基于浮架为刚体的假设，且水动力实验中网箱基本表现为刚体结构^[9-13]。在实际深远海海域，规格大的网箱在强浪流下会发生较大程度变形，深入了解浮架变形对网箱系统的可靠性设计至关重要。崔勇等^[14]基于有限元方法，对一种双层网底鲆鲽网箱位移与倾角进行了数值求解。Su等^[15]开发了一种结合数值模型和定位传感器数据的综合方法来实时监测网箱变形。刘海阳等^[16]采用数值模型分析了网箱部件在长期交变波浪荷载作用下导致的疲劳破坏问题。Chen等^[17]利用三维成像系统提出了一种估计网箱变形的方法。Tang等^[18]利用鱼群行为模型研究了鱼群对网箱流场和变形的影响。Dong等^[19]将浮架简化为圆环，基于曲梁理论分析了重力式网箱浮架系统的平面内外变形。Fredriksson等^[20]使用有限元建模方法预测了网箱浮架结构的临界载荷。Li等^[21]研究了半潜式网箱系统在波浪和水流中的动力响应。Huang等^[22]采用数值模型分析了均匀波流对多点系泊网箱系统的系泊索张力和网衣变形的影响。Kristiansen和Faltinsen^[23-24]提出了一种网衣黏性水动力荷载模型，研究了在水流和波浪作用下网箱的系泊荷载。

在海洋网箱养殖中，系泊方式主要包括单点系泊和多点系泊两种。多点系泊易于管理，但系泊成本较高^[25]；而单点系泊具有环境友好性，且系泊成本较低，但锚绳发生周期性疲劳破坏的机率相对更大^[26]。王绍敏等^[27]借鉴单点系泊系统，设计出一种可应用于往复流情况下的新型单锚腿系泊系统，并通过计算和实测证明了其有效性。梁家铭等^[28]设计了一种新型网箱浮架系统，与传统浮架进行了对比，并对4种常见的多点系泊系统进行了比较分析。我国南海海域已有部分养殖企业尝试应用单点系泊网箱，锚泊型式参照传统经验设计，缺乏科学依据。了解不同系泊方式深水网箱的受力变形特性，是解决网箱安全性能的前提。但国内外对于单多点系泊方式的水动力特性对比等研究相对较少。本研究针对一种实际应用的单点系泊高密度聚乙烯深水网箱，建立有限元数值计算模型，对原型网箱 (图1) 在波浪流作用下的浮架变形和锚绳受力进行了模拟研究，分析探讨了不同情况下系泊方式对网箱受力变形的影响，旨在为网箱浮架和系泊系统的优化设计提供理论依据和数据支撑。

图  1  三角形深水养殖网箱

Figure  1.  Scene of a triangular deep-water aquaculture cage

下载: 全尺寸图片 幻灯片

1.   数值模型

1.1   有限元模型

网箱系统主要由浮架、网衣和锚绳等主要部件组成，采用有限元方法将网箱系统离散为一系列线单元和集中质量点。集中质量点根据各部件实际情况分为三自由度质点 (图2) 和六自由度质点 (图3) 两种，均设于线单元两端。网箱的浮架系统离散为线单元和六自由度质点，其中六自由度质点起传递弯矩的作用。网衣和锚绳部件则采用线单元和三自由度质点模拟，考虑弯曲刚度很小，可设置为零。线单元所受的力和力矩包括浮力、波浪力、重力、剪力、张力和弯矩等，均集中分配到线单元两端的质点上，根据牛顿第二定律，建立如下质点运动方程：

图  2  三自由度质点

Figure  2.  3-D buoy

下载: 全尺寸图片 幻灯片

图  3  六自由度质点

Figure  3.  6-D buoy

下载: 全尺寸图片 幻灯片

式中M(p, a) 是惯性荷载，C(p, v) 是阻尼荷载，K(p) 为刚性荷载，F(p, v, t) 为外力荷载 (kN)；p为位置 (m)；v为速度 (m·s⁻¹)；a为加速度 (m·s⁻²)。

1.2   线单元模型理论

在波浪和水流共同作用下，网箱会产生运动和变形，因此高密度聚乙烯管制成的浮架所受的结构荷载包括因运动所产生的拉力、弯矩和剪力。根据网箱各部件实际情况，对于浮架的每个线单元，忽略扭转效应，对于网衣及锚绳等柔性部件，设定为只包含张力，抗弯刚度设为零。网箱系统不同部件的有效张力均可以表示为：

式中T_e表示有效张力 (kN)；T_w表示管壁张力 (kN)；P_i 和P₀分别为线单元内外部压力 (MPa)；A_i 和A₀分别为线单元内外横截面积 (m²)。式(3)右边第一项代表轴向刚度的作用，第二项是内外部压力的作用 (通过泊松比效应)，第三项代表轴向阻尼的作用。EA是轴向刚度 (kN)；ε为轴向应变，L为线单元的瞬时长度 (m)；λ为线单元的膨胀系数；L₀为线单元初始长度 (m)；μ为泊松比，对于高密度聚乙烯材料取值0.42； dL/dt是线单元长度增长率，C为阻尼系数，计算方式如下：

式中λ_a为目标轴向阻尼占比，取值0.4，选取这个值可以在不影响精度的情况下提供快速的模拟运行时间。C_crit为线单元的临界阻尼值，M表示线单元质量 (kg)。

对于网箱浮架结构所受的弯矩^[29-30]表示为：

式中EI为抗弯刚度 (kN·m²)；λ_b为目标弯曲阻尼；m为线单元质量 (kg)。

在求得各线单元端部的弯矩后，可以进一步计算线单元的剪力，各个线单元被看作是一个笔直的刚性杆，假设两端弯矩分别为M₁、M₂。由于刚性线单元的弯矩是沿着线单元线性变化的，故剪力为常数，即弯矩沿着线单元的变化率^[29-30]。

1.3   三自由度质点模型理论

三自由度质点仅具有升沉、纵移、纵摇3个自由度，由于不能发生旋转，主要用于模拟网衣连接点以及网衣底部沉子，可以模拟所受到的重力、浮力、阻力、附加质量的作用。

式中${F_{{\rm{D}}{\rm{rag}}}}$为阻力；${P_{\rm{W}}}$为浸水比例；ρ为水体密度；${{C}}_{{\rm{dx}}}$和${A_{\rm{x}}}$分别为x方向的阻力系数和阻力面积，${V_{\rm{r}}}$ 为流体相对于质点的速度矢量，$|{V_{\rm{r}}}|$为其绝对值, ${{V}}_{{\rm{rx}}}$为速度矢量在x轴上的分量；${F_{{\rm{Fluid}}}}$是流体加速度力，C_a为附加质量系数，A为流体加速度。

1.4   六自由度质点模型理论

六自由度质点具有3个方向的平动自由度和3个方向的旋转自由度，六自由度质点又分为集中质量点 (Lump buoy) 和圆柱式浮标 (Spar buoy)，集中质量点用于模拟浮架连接点。圆柱形浮标可以更好地模拟与水面的交互作用，用于模拟浮筒，通过计算瞬时位置和方向，可以确定传递力矩。2种六自由度质点的水动力荷载均使用莫里森方程计算：

式中$\rho $为海水密度，取1 025 kg·m⁻³；S为浮管单元投影面积 (m²)；V_r为水质点与线单元构件的相对速度 (m·s⁻¹)；$\forall _{\rm{w}}$为线单元构件排开水体的体积 (m³)；a_f为水质点加速度，a_r为相对加速度 (m·s⁻²)；C_d和C_a分别为速度力系数和附加质量系数。由于数值模拟采用的波浪理论为Ariy波，根据我国海港水文规范，经验系数C_d取值1.2，C_a取值1。

对于集中质量点，忽略水动力力矩，对于圆柱形浮标，受到的力矩为：

P_w为浮标的淹没比例，UDM表示的是平行或垂直于圆柱轴的单位阻尼力矩，流体假设无旋，$\varOmega $为浮标的相对角速度。

2.   模型验证

为了验证网箱数值模型的有效性及准确性，制作网箱整体模型比尺为1∶15的大比例尺模型，在中国水产科学研究院如东试验基地的波浪流水池中进行网箱水动力试验。水池长50 m、宽26 m、深1.2 m，配备了伺服电机驱动式造波机和造流系统，可产生试验要求的规则波、不规则波和稳态流，网箱模型参数见表1，网箱模型试验布置见图4。网箱浮架整体呈三角形状，内呈正六边形，网衣悬挂在正六边形的浮管上，其中对于锚绳长度和直径选取大比尺1∶15，使用与原材料相一致的聚乙烯绳。对于浮筒的几何比尺同样选取1∶15，根据Froude数相似律可以进一步推算出浮筒的浮力比尺为1∶15³，沉子重量同理。网衣的相似条件一般应遵循相似准则，但当模型比尺较大 (如λ=1∶15)时，不仅使网衣模型难以实现，而且会对网衣周围流体的流态发生较大改变，故参考Fredriksson等^[20]，在满足网衣主要水动力因素相似的条件下，对于网目大小和网线直径采用变比尺$ \lambda {\text{'}}$=1，对整体网衣三维空间尺寸则采用大比尺λ=1∶15。试验网衣材质与原型相同，均为PE。

表  1  网箱主要参数

Table  1.  Parameters of net cage

组件 Component	参数 Parameter	模型值 Model value	原型值 Prototype value
浮管 Floating pipe	直径 Diameter/mm	26.67	400
	壁厚 Wall thickness/mm	1.57	23.5
	形状 Shape	三角形	三角形
	长度 Length/m	1.33	20
网衣 Net	网目 Mesh size/mm	45	45
	目脚直径 Twine diameter/mm	2	2
	高度 Height/m	0.4	6
	材料 Material	PE	PE
锚绳 Mooring line	长度 Length/m	4	60
	直径 Diameter/mm	2.67	40
	材料 Material	PE	PE
浮筒 Buoy	高度 Height/m	0.13	2
	直径 Diameter/m	0.1	1.5
沉子 Sinker	质量 Mass/g	20×5.93	20×20000
	材料 Material	混凝土	混凝土

下载: 导出CSV 

| 显示表格

图  4  单点系泊网箱模型试验布置图

Figure  4.  Physical model configuration of SPM net cage

下载: 全尺寸图片 幻灯片

网箱系泊锚绳张力测量采用LA1型水下拉力传感器，该测力计设置于锚绳与网箱浮架的连接点A处位置，测力计量程为50 N，测试相对误差小于1%。采用非接触式运动姿态测量系统 (Untouched 6-D measurement system) 测量浮架中心点位置的运动响应情况。波浪测量采用LYL-Ⅲ型浪高仪，试验前后均对仪器标定校准，流速测量采用声学多普勒点式流速仪，测试相对误差小于0.5%。

试验波浪为规则波，波高设为16、20、24 cm，周期设为1.64、2.01 s；流速设为10.3、12.9 cm·s⁻¹。波浪和水流均沿x正方向传播。对于物理模型试验，试验过程中，一个工况结束后待水体静止之后再做下一个工况，各工况分别进行3次试验测量。图5为3种纯波条件下网箱最大锚绳张力、最大浮架升沉和最大浮架纵摇的计算结果与试验数据的对比情况，波浪周期为1.64 s。可以看出，3种不同波况条件下，网箱动力响应计算结果与模型试验结果较接近，相对误差均在10%以内。图6为4种波流组合条件下网箱锚绳力的时间历程曲线图。比较得知，最大锚绳张力的计算结果与试验结果相差较小，但两者曲线形状稍有不同，单个波浪周期内试验结果的曲线形状相比模拟结果的曲线稍尖些，这可能是由试验过程中锚绳力存在一定的脉冲效应造成的，使得锚绳力试验值在更短时间内达到峰值，但总体上两者曲线吻合较好，说明本研究建立的单点系泊网箱有限元模型可以模拟纯波、波流组合条件下的网箱动力响应情况。

图  5  纯波条件下网箱动力响应计算结果与试验值比较

Figure  5.  Comparison of dynamic responses between calculated and measured results of net cage under pure waves

下载: 全尺寸图片 幻灯片

图  6  波流条件下网箱锚绳力比较

Figure  6.  Comparison of mooring force for net cage under waves and current

下载: 全尺寸图片 幻灯片

3.   结果与分析

由于比尺缩放造成浮管刚度无法满足相似，采用物理模型试验测试波流条件下网箱浮架的变形难以实现，试验过程中的网箱浮架基本为刚体^[31]。因此，通过数值模拟方法研究网箱浮架的动力变形是切实可行的手段。深海区经常会有台风等恶劣天气，这种情况下波浪和洋流会同时出现；系泊系统将限制网箱的移动，以保证施工和养鱼的安全。在波浪和水流的共同作用下，锚缆的张力会在一定程度上反映网箱整体结构的受力。一旦系泊系统发生故障，不仅会造成鱼体逃逸和死亡，还会对网箱管理过程造成严重阻碍，给养殖人员带来巨大的经济损失，因此对系泊系统的安全性进行研究十分必要。基于本研究建立的有限元数值模型，在实际海况条件下，对原型尺寸的深水养殖网箱在不同系泊方式下进行水动力模拟分析，图7为单点系泊和多点系泊的布置方式，波流均沿x轴正向，对于多点系泊，3根锚绳呈120°角度分布。网箱各部件原型参数见表1。2种系泊方式的锚绳都是由浮筒按2∶1分为两段。设定的原型波浪流海况条件：波高为4、5、6 m，周期9 s，波浪为规则波；流速取0.5、1.0、1.5 m·s⁻¹；水深统一设为20 m。

图  7  不同系泊方式的网箱布置图

上部水平线为海平面，下部水平线为海底

Figure  7.  Layout of cage with different mooring systems

The upper level is sea level and the lower level is the seabed.

下载: 全尺寸图片 幻灯片

3.1   网箱系泊力比较

锚泊系统是网箱的重要组成部分，为了研究系泊方式对网箱系泊力的影响，对两种系泊方式在不同波流组合条件下进行了数值模拟。其中多点系泊仅考虑最危险工况，即单侧锚绳为迎浪流方向，考虑迎浪流侧系泊力。2种系泊方式下各工况对应的系泊力峰值计算结果见图8。可以看出，多点系泊方式系泊力峰值大于单点系泊方式，并且随着流速增加，这种差距更加明显。流速0.5 m·s⁻¹、波高4 m时，2种系泊方式系泊力峰值差值最小 (5.92%)，但随着流速和波高逐渐增加到更高工况后，两者差值为9.12%，差值增幅达到54.1%。这主要是由于背浪侧浮筒受到波流作用后，会通过锚绳向浮架施加一个向右的拉力，导致迎浪测的锚绳力比单点系泊情况下锚绳受力要大，这一点可从图9得到证明，观察多点系泊两侧锚绳力差值与单点系泊锚绳力历时曲线，可以看出除曲线谷值有所差异外，整体基本重合，说明增加的锚绳力主要来自于背浪侧锚绳。图9为波高6 m、流速1.5 m·s⁻¹时的锚绳历时曲线，由于背浪侧的两根锚绳呈对称布置，故仅考虑其中一根锚绳，可以看出当波流沿着x轴正向时，背浪侧系泊力峰值为8.06 kN，远小于迎浪测系泊力峰值(71.61 kN)，相差幅度达到788.46%。这是由于波流沿x轴正向时，背浪侧锚绳几乎不受力所致。

图  8  各工况下不同系泊方式的锚绳力峰值

Figure  8.  Peak value of mooring force under different sea conditions

下载: 全尺寸图片 幻灯片

图  9  锚绳力历时曲线 (H=6 m, v=1.0 m·s⁻¹)

Figure  9.  Time history of mooring force under wave height of 6 m and current velocity of 1.0 m·s⁻¹

下载: 全尺寸图片 幻灯片

3.2   浮架变形比较

由于高密度聚乙烯制成的浮架在大浪强流下会发生变形，将对网结构安全性能造成威胁，为了解2种系泊方式下浮架不同位置处的变形情况，且由于浮架为对称结构，仅选取沿y轴正向部分的多个点位进行分析，设置A—H点如图10所示，涉及的应力指标选取冯·米塞斯应力 (von Mises Stress)^[32]，当浮架线单元考虑扭转效应时，会计算出整个横截面上的应力峰值，应力计算表达式为：

图  10  浮架选定点位示意图

Figure  10.  Different points of floating collar

下载: 全尺寸图片 幻灯片

式中σ₁、σ₂、σ₃分别为3个方向的主应力 (MPa)。

表2列出了2种系泊方式下，大浪强流条件 (H=6 m, T=9 s, v=1.5 m·s⁻¹)浮架上A—H点所受的最大应力值。可以看出，2种系泊方式下的应力分布相对接近，单点系泊条件下浮架所受应力最大处为D点截面，达19.50 MPa，多点系泊时D点截面应力峰值为18.01 MPa，比单点系泊D点应力降低7.64%。其次为C点，应力峰值分别达14.99和13.88 MPa，相差7.4%。说明多点系泊能够一定程度上减少浮架的变形。但观察浮架D处达到应力峰值时，2种系泊方式下浮架受力变形都十分严重 (图11)，将对网箱的安全性能构成威胁，结合模拟计算结果接近HDPE材料的屈服强度(24 MPa)，若将网箱设置在更危险的海况下，需要采取合理措施对浮架D处进行加固处理。

表  2  浮架各点处变形

Table  2.  Deformation at different points of floating collar

浮架点编号 Point No.	应力 Stress/MPa
	单点系泊 SPM	多点系泊 MPM
A	11.07	11.77
B	3.44	2.93
C	14.99	13.88
D	19.50	18.01
E	2.34	2.68
F	2.11	3.03
G	3.16	2.25
H	7.37	3.71

下载: 导出CSV 

| 显示表格

图  11  D处达到应力峰值时浮架的变形情况

Figure  11.  Deformation of floating frame when stress peak reached at point D

下载: 全尺寸图片 幻灯片

3.3   网箱运动特性比较

3.3.1   浮架运动位置变化

为了进一步比较2种系泊方式下网箱在波流联合作用下的运动特性，首先进行了浮架的运动分析。取浮架上对称的A、D、E 3点进行运动分析，模拟工况为：波高6 m、周期9 s、流速1.5 m·s⁻¹。图12—图14分别为A、D、E点的x、y、z坐标历时曲线，由图可知，在波流作用下，系泊方式对浮架上A、B、C 3点的z坐标几乎无影响，说明改变系泊方式对浮架升沉影响不大；观察y坐标历时曲线，整体上看，主要差异在曲线峰谷值上，观察历时曲线稳定后，单点系泊条件下浮架D点的y坐标峰值为7.443 m，谷值为6.799 m，多点系泊条件下，y坐标峰值为7.416 m，谷值为6.876 m，峰谷差值降低了17.7%。说明多点系泊对浮架横摇有一定限制作用；对比A点x坐标，多点系泊时x坐标峰值为7.07 m，单点系泊时为5.627 m，相差25.64%。这是由于多点系泊条件下，背浪侧浮标在波流作用下产生向右的运动，从而对浮架产生向右的拉力，使得浮架在x轴上的位移增加，即多点系泊时浮架的纵移会更加明显。

图  12  A点坐标x、y、z历时曲线

Figure  12.  x, y and z duration curves of Point A

下载: 全尺寸图片 幻灯片

图  13  D点坐标x、y、z历时曲线

Figure  13.  x, y and z duration curves of Point D

下载: 全尺寸图片 幻灯片

图  14  E点坐标x、y、z历时曲线

Figure  14.  x, y and z duration curves of Point E

下载: 全尺寸图片 幻灯片

3.3.2   网衣运动位置比较

图15为网箱系统在流速1.5 m·s⁻¹、波高6 m、周期9 s条件下，一个完整的波周期内不同时刻网衣的运动状态，对比4个时刻网衣运动状态，变形情况几乎一致，说明改变系泊方式对此网箱的网衣体积几乎无影响。

图  15  一个波周期内网衣的xz面轮廓图

v=1.5 m·s⁻¹，H=6 m，T=9 s

Figure  15.  xz surface contour of net during a wave period

下载: 全尺寸图片 幻灯片

4.   结论

1)本文利用物理模型试验证明了网箱数值模型的有效性，分析了单点系泊与多点系泊系统网箱在系泊力、浮架变形以及网箱运动特性方面的区别。可以发现，对于不同的系泊方式网箱系统在波流作用下，不考虑背浪侧锚绳，多点系泊所受的系泊力峰值高于单点系泊，且随着流速增加，两者差异更为显著。

2)对于浮架的受力变形，两者浮架的应力峰值发生在相同位置，并且不同位置应力变化趋势一致，变形大致相同，但多点系泊浮架危险点处的应力变形低于单点系泊，两者应力峰值相差7.64%。

3)当模拟波高6 m、周期9 s、流速1.5 m·s⁻¹时，对浮架上选定的3个点位进行了x、y、z坐标的历时曲线分析，通过观察各点可以看出，对于此网箱系统，改变系泊方式对整个浮架的升沉运动几乎没有影响，但当波流沿x轴正向时，由于背浪侧浮筒在波流作用下产生向右的运动，进而导致多点系泊情况下浮架在x轴上的位移大于单点系泊，两者相差25.64%。进一步对比了网衣的运动状态，通过观察xz面在一个波周期内不同时刻网衣轮廓曲线，发现网衣的变化基本一致。说明改变系泊方式对于网衣系统的变形没有影响。