張 婧，華福永，周 游，施興華

（江蘇科技大學(xué)船舶與海洋工程學(xué)院，江蘇鎮(zhèn)江 212100）

0 引 言

近年來，隨著海洋工程領(lǐng)域有關(guān)技術(shù)和設(shè)備的不斷發(fā)展，海水養(yǎng)殖業(yè)不斷從近海走向遠海、深海，網(wǎng)漁具和網(wǎng)箱在海洋漁業(yè)生產(chǎn)中的應(yīng)用越來越廣泛。網(wǎng)箱上的網(wǎng)衣作為其主要構(gòu)件，對整體設(shè)備的受力起主導(dǎo)作用。與其他海洋結(jié)構(gòu)物相同，網(wǎng)衣及由網(wǎng)衣構(gòu)成的漁具和漁業(yè)設(shè)施在海洋中會承受風、浪和流等環(huán)境因素的影響。網(wǎng)衣的水動力研究［1-3］涉及多個學(xué)科，包括數(shù)學(xué)、物理學(xué)和流體力學(xué)等，研究手段非常復(fù)雜。

目前有關(guān)網(wǎng)衣類型對流場特性和水動力性能的影響的研究尚未十分成熟。在國外方面：PATURSSON［4］提出采用多孔介質(zhì)模型模擬平面網(wǎng)衣，研究了水流經(jīng)過網(wǎng)衣作用之后，其速度衰減的規(guī)律；JOHANSSON等［5］測量了在挪威西海岸區(qū)域和在峽灣內(nèi)布置的組合式網(wǎng)箱周圍的流速，分析了網(wǎng)箱的流速衰減情況，同時對網(wǎng)箱周圍的鹽度、溫度和溶解氧進行了測量；WERGELAND［6］研究了網(wǎng)箱配重底圈的渦激振動問題，分析了組合式網(wǎng)箱在波浪作用下的錨繩受力特性和在波流聯(lián)合作用下的體積衰減問題；PASCAL 等［7］采用聲學(xué)多普勒流速剖面儀和聲學(xué)多普勒流速儀（Acoustic Doppler Velocimeter，ADV）測量了水流衰減情況，研究了水流經(jīng)過圓形網(wǎng)箱的作用之后，其速度衰減的規(guī)律。在國內(nèi)方面：劉超等［8］設(shè)計了4 種不同網(wǎng)目形狀和不同網(wǎng)目大小的網(wǎng)箱模型，采用ADV和粒子圖像測速（Particle Image Velocimetry，PIV）技術(shù)對網(wǎng)箱系統(tǒng)內(nèi)外的流速和流場變化進行了測定，研究了圓形重力式網(wǎng)箱在動態(tài)流場中的水動力學(xué)性能，同時獲得了網(wǎng)箱系統(tǒng)內(nèi)部與周圍的流場分布及流速變化規(guī)律；林為付［9］通過開展模擬試驗研究了與網(wǎng)箱動態(tài)響應(yīng)有關(guān)的參數(shù)，分析了圓形養(yǎng)殖網(wǎng)箱結(jié)構(gòu)及流場中網(wǎng)箱的動態(tài)響應(yīng)情況；李鵬等［10］基于計算流體力學(xué)（Computational Fluid Dynamics，CFD）方法對小尺度平面無結(jié)金屬網(wǎng)衣進行了數(shù)值模擬試驗，研究了其在均勻流作用下的水動力性能和網(wǎng)衣周圍的流場分布，得出了阻力系數(shù)和速度場分布隨雷諾數(shù)、網(wǎng)衣的攻角、網(wǎng)衣密實度和網(wǎng)眼形狀的變化趨勢。

在水流作用下，網(wǎng)衣受力會使網(wǎng)箱發(fā)生變形，一旦水流速度過大，不僅會導(dǎo)致網(wǎng)箱的容積大幅度減小，而且會影響網(wǎng)箱內(nèi)外水體的交換，降低養(yǎng)殖生物的存活率?；谏鲜銮闆r，本文通過STARCCM ＋軟件進行數(shù)值仿真，對比分析雙平面網(wǎng)衣與圓形網(wǎng)衣在流場中不同流速下的拖曳力和阻力系數(shù)的變化規(guī)律，為海洋網(wǎng)箱設(shè)計和流場特性分析提供參考。

1 數(shù)值模型

1.1 網(wǎng)衣數(shù)值模型

基于三維建模軟件SolidWorks建立雙平面網(wǎng)衣和圓形網(wǎng)衣2 種網(wǎng)衣模型，具體參數(shù)見表1。為著重研究網(wǎng)衣在遮蔽效應(yīng)下的受力情況，將方形網(wǎng)箱上的網(wǎng)衣簡化為同一方向上的前后2 個受力平面，將圓形網(wǎng)箱外網(wǎng)衣簡化為水平流域上的受力面。雙平面網(wǎng)衣模型參照方形網(wǎng)箱模型及其不相鄰的2 個平面，雙網(wǎng)衣之間的距離為0.8 m。圓形網(wǎng)衣模型參照圓形網(wǎng)箱的形狀，為確保其面積與平面網(wǎng)衣的投影面積相同，取圓形網(wǎng)衣的直徑為0.8 m，高為0.4 m，其他條件與平面網(wǎng)衣基本保持一致。

表1 網(wǎng)衣參數(shù)

1.2 網(wǎng)衣受力計算

網(wǎng)衣是多孔狀結(jié)構(gòu)，一般采用水池試驗或數(shù)值計算的方法對其水動力性能進行研究；計算網(wǎng)衣流體載荷的方法主要有多孔介質(zhì)法、集中質(zhì)量法、屏模型和莫里森公式等。下面采用莫里森公式對2 種網(wǎng)衣進行數(shù)值計算，包括網(wǎng)衣受力面改進計算和水動力系數(shù)選取，在此基礎(chǔ)上建立流場域，設(shè)置數(shù)值模擬的邊界條件。

1.2.1 網(wǎng)衣受力相似條件

對于網(wǎng)衣而言，控制受力的動力因素主要是水流阻力、網(wǎng)衣在水中的重量和浮力，其中網(wǎng)衣所受水流阻力與網(wǎng)衣的垂直投影面積有關(guān)，這樣從相似準則的角度來說，只要滿足主要動力因素的相似要求即可。一般而言，網(wǎng)衣受力相似條件應(yīng)遵循幾何相似、水阻力相似、重量和浮力相似。本文在已有研究的基礎(chǔ)上，提出以下計算改進方法：在保持網(wǎng)衣總尺度幾何相似的前提下，確保網(wǎng)衣受力總面積幾何相似，增大網(wǎng)線直徑和網(wǎng)目尺寸，減小模型網(wǎng)衣對水流流態(tài)的改變，保證水動力相似并使網(wǎng)衣模型可行。

在計算網(wǎng)衣受力時，可采用流體動力的一般公式（Morison方程），即

式（1）中：F為網(wǎng)衣的水阻力；A為網(wǎng)衣與水流垂直時的網(wǎng)線投影面積；CD為水阻力系數(shù)；ux為流速；ρ為流體密度。當考慮有結(jié)節(jié)網(wǎng)衣時，結(jié)節(jié)直徑按3.14 乘以目腳直徑計算，單位縮結(jié)面積的網(wǎng)線投影面積為

若忽略式（2）右端的第二項，并取整個網(wǎng)衣的縮結(jié)面積為S，則式（2）可近似表示為

日本學(xué)者田內(nèi)［11］曾考慮各種沖角，并對網(wǎng)衣進行阻力測試，得到經(jīng)驗公式

式（2）～式（4）中：a為目腳；θ為沖角；S為網(wǎng)衣的大尺度面積。當θ ＝90°時，式（4）可簡化為

1.2.2 網(wǎng)衣的水動力系數(shù)選取

增大網(wǎng)線直徑和網(wǎng)目尺寸會導(dǎo)致網(wǎng)衣體積增大，產(chǎn)生的浮力遠大于改進之前的浮力，因此為協(xié)調(diào)重力與浮力，需調(diào)整網(wǎng)衣的水動力系數(shù)。在純水流情況下，網(wǎng)衣的水動力系數(shù)CD與雷諾數(shù)Ren有關(guān)，將水阻力系數(shù)與雷諾數(shù)的關(guān)系總結(jié)成函數(shù)關(guān)系，有

式（6）～式（8）中：CDn為速度力的切向速度力系數(shù)；CDt為速度力的法向速度力系數(shù)；μ 為液體的運動黏性系數(shù)；ρ為流體密度；vRn為水流相對目腳的流速；D為網(wǎng)目目腳直徑。液體的運動黏性系數(shù)與溫度有關(guān)，這里取溫度為20 ℃時的值，即1.35 ×10－6m2／s。為方便計算，將速度力系數(shù)設(shè)定為常數(shù)，由于CDn在1.1 附近小范圍波動，CDt在0.14 附近小范圍波動，因此分別取CDn和CDt的值為1.1 和0.14。

1.3 控制方程

流體流動連續(xù)性方程［11］的微分形式為

對于密度為常數(shù)的不可壓縮均質(zhì)流體，式（10）可進一步簡化為

式（10）和式（11）中：ρ為流體的密度；t為時間；u、v和w分別為速度矢量在x方向、y方向和z方向的分量。

2 網(wǎng)格收斂性驗證及網(wǎng)格劃分

2.1 網(wǎng)格收斂性驗證

設(shè)定網(wǎng)衣數(shù)值模擬的邊界條件，對全流域網(wǎng)格的基礎(chǔ)尺寸進行調(diào)整，劃分出數(shù)量不同的3 套網(wǎng)格，由此對網(wǎng)格的收斂性進行驗證。3 套網(wǎng)格的基礎(chǔ)尺寸設(shè)置及網(wǎng)格總數(shù)見表2。

表2 3 套網(wǎng)格的基礎(chǔ)尺寸設(shè)置及網(wǎng)格總數(shù)

3 套網(wǎng)格的截圖見圖1。從圖1 中可看出：隨著網(wǎng)格數(shù)量的增加，網(wǎng)衣表面和計算域的網(wǎng)格劃分越來越均勻，同時兩者交界處的過渡越來越自然；隨著網(wǎng)格數(shù)量的增加，計算的負擔越來越大。因此，為選擇合適的網(wǎng)格尺寸，以進口速度0.5 m／s為例對3 套網(wǎng)格進行后網(wǎng)衣阻力計算，分別對比達到穩(wěn)定所需迭代步數(shù)和阻力值。

圖1 3套網(wǎng)格的截圖

圖2 為3 套網(wǎng)格的迭代步數(shù)與阻力的關(guān)系曲線對比。從迭代曲線上看，3 套網(wǎng)格的迭代步數(shù)到400 步時就已趨于平穩(wěn)；從運算結(jié)果上看，3 套網(wǎng)格的受力值相差不到0.5 N，相對誤差較小，在6%以內(nèi)。綜上所述，考慮到網(wǎng)格過多會影響計算機的運算速度，后續(xù)水動力數(shù)值模擬分析均采用基礎(chǔ)尺寸為250 mm 的網(wǎng)格進行有關(guān)計算。

圖2 3套網(wǎng)格的迭代步數(shù)與阻力的關(guān)系曲線對比

2.2 網(wǎng)格劃分

網(wǎng)格是幾何域的離散表示，也是計算模型前處理的重要組成部分，因此本文建立計算域，并對其進行網(wǎng)格劃分。將網(wǎng)衣置于尺寸為5.0 m×2.0 m×0.6 m 的計算域內(nèi)，網(wǎng)格生成器為表面重構(gòu)和切割體網(wǎng)格單元生成器，網(wǎng)格基礎(chǔ)尺寸為0.25 m，最小表面尺寸為基礎(chǔ)尺寸的10%（即0.025 m）。在計算域內(nèi)布置尺寸為5.0 m×2.0 m×0.2 m的自由液面，為對自由液面網(wǎng)格進行加密，設(shè)置網(wǎng)格各向異性，X軸和Y軸方向的網(wǎng)格大小相當于基礎(chǔ)尺寸的20%（為0.050 0 m），Z軸方向的網(wǎng)格大小相當于基礎(chǔ)尺寸的5%（為0.012 5 m）。

3 網(wǎng)衣三維數(shù)值模擬

3.1 流場模擬

選取入射水流速度分別為0.5 m／s、1.0 m／s、1.5 m／s和2.0 m／s，對網(wǎng)衣的整體流場進行數(shù)值模擬，入射角為90°，模擬時間為20 s，隱式不定常時間為0.008 s，最大內(nèi)部迭代次數(shù)為10 次。著重分析2 種網(wǎng)衣在水平流場中的變化。圖3 和圖4 為2 m／s流速下2 種網(wǎng)衣的水平流場分布。

圖3 2 m／s流速下雙平面網(wǎng)衣的水平流場分布

圖4 2 m／s流速下圓形網(wǎng)衣的水平流場分布

從圖3 和圖4 中可看出：由于網(wǎng)衣的遮蔽效應(yīng)，水流經(jīng)過2 種網(wǎng)衣時均會產(chǎn)生阻流，在網(wǎng)衣之間均存在明顯的流速衰減；當水流穿過2 層網(wǎng)衣之后，在后網(wǎng)衣之后的區(qū)域流速衰減幅度更大，在圖3 中靠近后網(wǎng)衣的區(qū)域以及圖4 中靠近兩側(cè)網(wǎng)衣的區(qū)域存在一段很短的壓力真空帶；隨著與網(wǎng)衣距離的增大，流速衰減區(qū)域的寬度有逐漸減小的趨勢；衰減區(qū)域兩側(cè)及底部的水流速度，隨著與網(wǎng)衣距離的增加而增大；由于網(wǎng)衣的開孔率較大，當水流穿過網(wǎng)衣時，其流動方向并沒有發(fā)生明顯的變化。

3.2 網(wǎng)衣水平拖曳力分析

進一步對4 種流速場下2 種網(wǎng)衣的受力變化進行對比分析。4 種流速下雙平面網(wǎng)衣和圓形網(wǎng)衣的受力曲線分別見圖5 和圖6。

圖5 4種流速下雙平面網(wǎng)衣的受力曲線

圖6 圓形網(wǎng)衣受力曲線

由圖5 可知，在流場中，雙平面網(wǎng)衣的水平拖曳力在初始階段存在一個極大值，從極大值開始受力迅速減小，隨后迅速增大，最后達到穩(wěn)態(tài)，穩(wěn)態(tài)受力結(jié)果見表3。在同一流速下，前網(wǎng)衣達到穩(wěn)態(tài)的時間比后網(wǎng)衣早，前網(wǎng)衣受到的水平拖曳力比后網(wǎng)衣大。對比穩(wěn)態(tài)受力結(jié)果可知，隨著流速的增加，網(wǎng)衣的水平拖曳力增大，前網(wǎng)衣與后網(wǎng)衣的受力結(jié)果相差越來越大。圖7 為2 種網(wǎng)衣在穩(wěn)態(tài)下的水平拖曳力對比。由圖7 可知，圓形網(wǎng)衣在流場中的受力變化與雙平面網(wǎng)衣相似，當流速增加時，網(wǎng)衣的水平拖曳力增大。

圖7 2種網(wǎng)衣在穩(wěn)態(tài)下的水平拖曳力對比

表3 雙平面網(wǎng)衣穩(wěn)態(tài)受力結(jié)果

由此可知，圓形網(wǎng)衣的受力始終比雙平面網(wǎng)衣的大，兩者在流速為0.5 m／s時相差最小，為5.3 N，在流速為2 m／s時相差最大，為48.6 N。圓形網(wǎng)衣與雙平面網(wǎng)衣的水平拖曳力有相同的變化趨勢，均隨著流速的增加而增大。

3.3 阻力系數(shù)分析

根據(jù)遮蔽效應(yīng)下的流速衰減和網(wǎng)衣受力情況，計算不同流速下網(wǎng)衣的阻力系數(shù)。4 種流速下雙平面網(wǎng)衣和圓形網(wǎng)衣的阻力系數(shù)隨時間的變化曲線分別見圖8 和圖9。

圖8 雙平面網(wǎng)衣的阻力系數(shù)隨時間的變化曲線

圖9 圓形網(wǎng)衣的阻力系數(shù)隨時間的變化曲線

由圖8 可知：在流場中，雙平面網(wǎng)衣的阻力系數(shù)在初始階段迅速減小，隨后迅速增加，最后達到穩(wěn)態(tài)；在同一流速下，前網(wǎng)衣達到穩(wěn)態(tài)的時間比后網(wǎng)衣早。整體上看，前網(wǎng)衣的阻力系數(shù)比后網(wǎng)衣的略大；當流速為1.0 m／s時，前網(wǎng)衣的阻力系數(shù)比后網(wǎng)的略小，但因該流速下兩者的結(jié)果基本一致，對網(wǎng)衣的影響較小，可忽略不計。對比穩(wěn)定狀態(tài)下前后網(wǎng)衣的阻力系數(shù)值發(fā)現(xiàn)，前網(wǎng)衣與后網(wǎng)衣的阻力系數(shù)相近，當流速為1.5 m／s時，兩者的阻力系數(shù)僅相差0.004，當流速增加時，雙平面網(wǎng)衣的阻力系數(shù)反而減小。

圖10 為2 種網(wǎng)衣的阻力系數(shù)隨流速的變化曲線對比。由圖10 可知：穩(wěn)態(tài)下圓形網(wǎng)衣與雙平面網(wǎng)衣的阻力系數(shù)相差不大，且隨著流速的增加，兩者的數(shù)值越來越接近；當流速為0.5 m／s時，2 種網(wǎng)衣的阻力系數(shù)相差最大，為0.066；當流速為2.0 m／s時，2 種網(wǎng)衣的阻力系數(shù)僅相差0.008。與雙平面網(wǎng)衣一樣，當流速增加時，圓形網(wǎng)衣的阻力系數(shù)反而減小。

圖10 2種網(wǎng)衣的阻力系數(shù)隨流速的變化曲線對比

4 結(jié) 語

本文運用STARCCM ＋軟件對2 種網(wǎng)衣進行數(shù)值仿真模擬，建立合適的流場域驗證網(wǎng)格的收斂性。通過分析水流經(jīng)過網(wǎng)衣之后的流速衰減情況，計算網(wǎng)衣的受力和阻力系數(shù)，主要得到以下結(jié)論：

1）由于網(wǎng)衣的遮蔽效應(yīng)，在網(wǎng)內(nèi)和網(wǎng)后存在明顯的流速衰減，在靠近網(wǎng)衣的區(qū)域可清晰地看見一段很短的真空帶，隨后流速增加，逐漸達到穩(wěn)態(tài)。在水平流場中，隨著與網(wǎng)衣距離的增大，在衰減區(qū)域兩側(cè)水流的速度相應(yīng)增加。

2）網(wǎng)衣網(wǎng)前流速均大于網(wǎng)內(nèi)和網(wǎng)后流速，網(wǎng)內(nèi)流速大于網(wǎng)后流速，各流速下網(wǎng)內(nèi)的流速衰減系數(shù)均大于網(wǎng)后的流速衰減系數(shù)，且都在一定范圍內(nèi)。

3）在流場中，網(wǎng)衣的水平拖曳力和阻力系數(shù)在初始階段存在一個極大值，從極大值開始受力迅速減小，隨后迅速增加，最后達到穩(wěn)態(tài)。當流速增加時，網(wǎng)衣的水平拖曳力增大，阻力系數(shù)減小。

4）圓形網(wǎng)衣的受力始終比雙平面網(wǎng)衣的大，隨著流速的增加，兩者具有相同的變化趨勢；穩(wěn)態(tài)下圓形網(wǎng)衣與雙平面網(wǎng)衣的阻力系數(shù)相差不大，且隨著流速的增加，兩者的數(shù)值趨近。