潘 昀，劉 燦，楊 熙，施圣鈐，沈良朵，桂福坤

(1. 浙江海洋大學(xué) 港航與交通運(yùn)輸工程學(xué)院,浙江 舟山 316022; 2. 浙江海洋大學(xué) 國家海洋設(shè)施養(yǎng)殖工程技術(shù)研究中心,浙江 舟山 316022)

南方海洋牧場環(huán)境因子不同于北方海域，并且主要以魚、蝦、蟹、頭足類等高級消費(fèi)者為增養(yǎng)殖對象，因此海洋牧場建設(shè)的模式應(yīng)以浮式魚礁、養(yǎng)殖筏式、養(yǎng)殖網(wǎng)箱等中上層人工生境構(gòu)建為基礎(chǔ)[1]，增殖天然餌料、養(yǎng)殖藻場和創(chuàng)造魚類保育場，為魚類產(chǎn)卵和幼體庇護(hù)提供良好的棲息環(huán)境。農(nóng)業(yè)部印發(fā)的《國家級海洋牧場示范區(qū)建設(shè)規(guī)劃(2017—2025)》中提出：截止到2025年，規(guī)劃在東海區(qū)共建設(shè)20個(gè)國家級海洋牧場示范區(qū)(包括2015—2016年已建情況)，形成示范海域面積500多平方千米，其中建設(shè)人工魚礁區(qū)面積160平方千米，投放人工魚礁500多萬空立方米?？梢?，國家高度重視和大力支持東海沿海地區(qū)海洋牧場建設(shè)，而如何選擇浮魚礁類型、確定投放布局、保證結(jié)構(gòu)安全和材料安全、保證聚魚和庇魚效果是科學(xué)創(chuàng)造魚類產(chǎn)卵場和保育場的重要前提。

日本是最早研究和使用中上層浮魚礁的國家，指出浮魚礁一般由浮力系統(tǒng)、礁體系統(tǒng)、錨泊系統(tǒng)三部分組成[2-3]，日本和韓國近年來也設(shè)計(jì)和投放了用于近海的典型中上層浮魚礁結(jié)構(gòu)型式[4-5]。浮力系統(tǒng)構(gòu)件主要為浮管和浮球組合而成，提供浮魚礁的浮力同時(shí)也是波浪荷載的主要來源；礁體部分一般由剛性鐵或者塑料框架覆蓋網(wǎng)衣、塑料薄膜、帆布構(gòu)成，網(wǎng)衣的透水性和柔性較大，且宜附著藻類生物為魚群營造產(chǎn)卵場所和索餌區(qū)，更受青睞，也是荷載的來源；不同的錨泊方式對浮魚礁的運(yùn)動(dòng)特性影響較大，需要搭配最佳的浮力和配重才能更好的發(fā)揮對浮體的系留作用，各系統(tǒng)之間及構(gòu)件的協(xié)調(diào)、搭配非常關(guān)鍵。中上層浮魚礁投放于近海，屬于海岸工程范疇，因此，進(jìn)一步開展中上層浮魚礁結(jié)構(gòu)安全、性能穩(wěn)定、系統(tǒng)構(gòu)成合理搭配的研究以及其抗浪機(jī)制的研究，對設(shè)計(jì)和投放人工浮魚礁具有重要的意義。

近年來用于海洋漁業(yè)增殖的浮魚礁也出現(xiàn)了新的結(jié)構(gòu)型式[6]，并進(jìn)行少數(shù)單體投放和聚魚效果評估，但針對海洋水動(dòng)力學(xué)特性方面的安全和穩(wěn)定性能以及浮魚礁群布設(shè)方式等諸多工程基礎(chǔ)問題和技術(shù)均少有報(bào)道和論文發(fā)表。浮魚礁的浮力系統(tǒng)主要由浮球或浮管排列組合而成，魚礁主體由綱繩和網(wǎng)衣編織而成，結(jié)構(gòu)型式與養(yǎng)殖網(wǎng)箱非常類似。Hou等[7]、Huang等[8]、Bai等[9]、Zhao等[10]、Gui等[11]、Moe-F?re等[12]、Strans等[13]長期致力于網(wǎng)箱結(jié)構(gòu)型式、網(wǎng)箱浮架運(yùn)動(dòng)、錨系方式、網(wǎng)箱容積變化的研究，才解決了網(wǎng)箱柔性結(jié)構(gòu)在波浪和潮流作用下體積縮小及結(jié)構(gòu)安全問題，確保網(wǎng)箱在不同海洋環(huán)境下滿足設(shè)計(jì)要求。這些均為波浪作用下的浮魚礁動(dòng)力響應(yīng)研究提供了保障，與網(wǎng)箱結(jié)構(gòu)形式不同的是中上層浮魚礁必須保持一定的高度和姿態(tài)才能作用于周圍水體的魚群[4]，垂向尺寸往往比水平尺寸大5～10倍，并且中上層浮魚礁多使用單錨系泊[4, 14-15](可群化投放)。數(shù)值建模方法與之不同的是復(fù)雜浮體結(jié)構(gòu)轉(zhuǎn)動(dòng)慣量及轉(zhuǎn)動(dòng)變換布萊恩角的求算，以及中上層浮魚礁各構(gòu)成組件連接位置處點(diǎn)或單元的相互耦合受力和運(yùn)動(dòng)處理方法，這是研究的重點(diǎn)。

1 數(shù)值模型建立

中上層浮魚礁各組成要素之間通過熱熔和系縛方式進(jìn)行連接，因此整體浮魚礁的運(yùn)動(dòng)數(shù)值模型建立的關(guān)鍵就是連接位置處點(diǎn)或單元的處理。首先介紹圓管浮架、網(wǎng)衣、錨繩等組成要素單獨(dú)波浪作用下受力和運(yùn)動(dòng)的數(shù)值計(jì)算，而后重點(diǎn)介紹三者之間相互耦合受力和運(yùn)動(dòng)的處理方法。

波浪作用下浮架、網(wǎng)衣及錨繩運(yùn)動(dòng)數(shù)值模擬方法分剛性和柔性兩種情況，浮架為剛性結(jié)構(gòu)，運(yùn)動(dòng)包括平動(dòng)和轉(zhuǎn)動(dòng)，而網(wǎng)衣和錨繩為柔性結(jié)構(gòu)，運(yùn)動(dòng)形式表現(xiàn)為拉伸和懸浮運(yùn)動(dòng)。數(shù)值模型主要采用莫里森方程計(jì)算劃分單元和集中質(zhì)量點(diǎn)的受力，進(jìn)而獲得整體結(jié)構(gòu)的合力和繞質(zhì)心的力矩，而后離散時(shí)間獲得結(jié)構(gòu)的運(yùn)動(dòng)狀態(tài)和受力變化。

1.1 圓管浮架受力和運(yùn)動(dòng)

浮架采用HDPE管道熱熔焊接，在小尺度范圍內(nèi)可認(rèn)為管道為剛性，能夠維持設(shè)計(jì)形狀。由于波浪沿水深方向水質(zhì)點(diǎn)速度和水質(zhì)點(diǎn)加速度逐漸減小，因此以浮體的質(zhì)心為中心，上部所受水平波浪力大于底部所受水平波浪力，這將導(dǎo)致浮魚礁礁體發(fā)生轉(zhuǎn)動(dòng)。對浮架進(jìn)行有限單元?jiǎng)澐趾笤儆?jì)算單元的波浪力，能夠準(zhǔn)確刻畫這一現(xiàn)象。

1.1.1 浮架受力和單元?jiǎng)澐?/p>

圓管浮架在波浪作用下的主要受力包括重力、浮力、波浪力以及網(wǎng)衣和錨繩的拉力。當(dāng)圓管的直徑與波長的比值小于0.2時(shí)，波浪力可使用莫里森方程計(jì)算。至于運(yùn)動(dòng)的圓柱桿件波浪力(包括拖曳力和慣性力)的計(jì)算以及水動(dòng)力系數(shù)的選取在文獻(xiàn)[16]中有較詳細(xì)的介紹。當(dāng)浮架發(fā)生轉(zhuǎn)動(dòng)后，波浪水質(zhì)點(diǎn)速度方向、水質(zhì)點(diǎn)加速度方向與圓管桿件之間的夾角發(fā)生變化，需要在圓管單元上建立局部坐標(biāo)系，具體可參見文獻(xiàn)[16]。

與文獻(xiàn)[16]不同的是，此次建立的是三維浮架運(yùn)動(dòng)模型并且其周圍纏繞網(wǎng)衣，如圖1和圖2所示。以浮架的正視面ABCD為例，為了保證網(wǎng)衣網(wǎng)目恰好完整地系縛在圓管上，水平方向AB和CD及AD和BC圓管的單元?jiǎng)澐中枰獫M足網(wǎng)目的尺寸大小。陳天華等[17-18]常采用水平縮結(jié)系數(shù)0.66(網(wǎng)目水平對角線長度與目腳長度2a=8 cm的比值)和垂直縮結(jié)系數(shù)0.75(網(wǎng)目垂直對角線長度與目腳長度2a=8 cm的比值)的網(wǎng)目形狀。如果網(wǎng)目群化技術(shù)[17-18]采用8×8，那么網(wǎng)目的水平長度為0.422 4 m和垂直長度0.480 0 m。圖1給出了框架浮體一個(gè)面的劃分單元方式，網(wǎng)目系縛在圓管的奇數(shù)單元上，因此水平方向單元的長度為0.211 2 m，垂直方向單元的長度為0.240 0 m。圖1中，每個(gè)方格表示一半單元長度的物理量(單元質(zhì)量、波浪場位置處的水質(zhì)點(diǎn)速度和水質(zhì)點(diǎn)加速度、計(jì)算力臂的位置等)集中點(diǎn)。因此在三維模型中，除四個(gè)端點(diǎn)(A、B、C、D)處相鄰三個(gè)單元結(jié)點(diǎn)，其余單元結(jié)點(diǎn)僅左右相鄰兩個(gè)其他單元結(jié)點(diǎn)。

圖1 浮架單元?jiǎng)澐旨捌涫芰κ疽釬ig. 1 Division and force analysis of the floating structure

圖2 浮架與網(wǎng)衣計(jì)算方法Fig. 2 Calculation method of the floating structure and netting

模型中，整體坐標(biāo)原點(diǎn)建立在浮架的中軸線上，z=0與靜水面重合，x、y和z軸的正方向如圖2(a)所示。給定浮架的初始位置坐標(biāo)分三種情況：一種是與z軸平行的圖2(a)中圓形標(biāo)號1～4號圓管，計(jì)算程序編寫的數(shù)組編號方向?yàn)樽陨隙?，如P1(ibar,4,3)，其中ibar表示劃分的ibar個(gè)單元點(diǎn)，4表示此類圓管的編號，3表示三維整體坐標(biāo)x、y和z值；另一種是與x軸平行的圖2(a)中方形標(biāo)號1～4號圓管，數(shù)組編號的方向?yàn)樽郧岸?，如P2(jbar,4,3)，其中jbar表示把此類圓管劃分成了jbar個(gè)單元點(diǎn)；第三種是與y軸平行的圖2(a)中三角形標(biāo)號1～4號圓管，數(shù)組編號的方向?yàn)樽宰蠖遥鏟3(jbar,4,3)，數(shù)字意義與前述一致。

1.1.2 浮架運(yùn)動(dòng)

浮魚礁的浮架周圍系縛纏繞網(wǎng)衣，在4根輔助錨和1根主錨的牽引下，浮架的運(yùn)動(dòng)包括平動(dòng)和轉(zhuǎn)動(dòng)。根據(jù)剛體運(yùn)動(dòng)的速度和加速度合成原理，對于剛體在固定坐標(biāo)系(絕對坐標(biāo)系)中的運(yùn)動(dòng)速度和加速度實(shí)際上就是牽連速度和牽連加速度，包括平動(dòng)和轉(zhuǎn)動(dòng)兩部分。

1) 計(jì)算平動(dòng)

先不考慮浮架轉(zhuǎn)動(dòng)，在整體坐標(biāo)系下將浮架所有單元的受力進(jìn)行累加，如圖1主要包括波浪力、網(wǎng)拉力、浮力、重力，而后根據(jù)牛頓第二定律獲得浮架的加速度，如式(1)：

(1)

式中：M表示浮架的總質(zhì)量，ΔM表示浮架各單元慣性質(zhì)量和(某一單元的體積為Vvolume，取拖曳力系數(shù)為CM，則慣性質(zhì)量為Δm=ρwater(CM-1)·Vvolume)；aframe1表示框架浮體的平動(dòng)加速度；T表示錨繩對框架的拉力；FD、FI分別表示浮架的合成拖曳力和慣性力；G表示浮架的重力；Ffloat表示浮架的浮力；Vframe表示浮架的平動(dòng)速度，V0為框架浮體的初始速度或上一時(shí)刻的速度。

2) 浮架轉(zhuǎn)動(dòng)慣量公式推導(dǎo)

正確計(jì)算浮架的轉(zhuǎn)動(dòng)慣量，是計(jì)算其轉(zhuǎn)動(dòng)的前提。如圖2(a)坐標(biāo)系，相對于整體坐標(biāo)系原點(diǎn)而言，以剛體的質(zhì)心作為動(dòng)坐標(biāo)的原點(diǎn)可進(jìn)行轉(zhuǎn)動(dòng)計(jì)算。根據(jù)平行軸定理，把浮架的12根圓管的轉(zhuǎn)動(dòng)慣量分x、y、z三個(gè)方向平移至質(zhì)心處。如圖2(a)所示，構(gòu)成浮架圓管可分為三類，即圓形標(biāo)號1、方形標(biāo)號1、三角形標(biāo)號1。假設(shè)圓管單位長度的質(zhì)量為m，圓標(biāo)號1劃分為ibar個(gè)單元結(jié)點(diǎn)，每個(gè)單元的長度為Lvertical，方形標(biāo)號1和三角形標(biāo)號1均劃分為jbar個(gè)單元，每個(gè)單元的長度為Llever，那么三類圓管的轉(zhuǎn)動(dòng)慣量平移至浮架質(zhì)心處的計(jì)算表達(dá)式為式(2)～(3)。

圓形標(biāo)號1圓管的x、y、z三方向轉(zhuǎn)動(dòng)慣量平移至質(zhì)心處：

(2)

方形標(biāo)號1圓管的x、y、z三方向轉(zhuǎn)動(dòng)慣量平移至質(zhì)心處：

(3)

三角形標(biāo)號1圓管的x、y、z三方向轉(zhuǎn)動(dòng)慣量平移至質(zhì)心處：Ix3=Iy2，Iy3=Ix2，Iz3=Iy3。因此，浮架總的x、y、z三方向轉(zhuǎn)動(dòng)慣量為：Ix=4(Ix1+Ix2+Ix3)，Iy=4(Iy1+Iy2+Iy3)，Iz=4(Iz1+Iz2+Iz3)。式(2)和式(3)中的rbar為浮架所用圓管的半徑，所有物理量均采用國際單位制計(jì)算。

3) 計(jì)算轉(zhuǎn)動(dòng)

剛體的轉(zhuǎn)動(dòng)是相對于整體坐標(biāo)系原點(diǎn)而言，一般常取整體坐標(biāo)系剛體的質(zhì)心作為動(dòng)坐標(biāo)的原點(diǎn)進(jìn)行轉(zhuǎn)動(dòng)計(jì)算。具體計(jì)算轉(zhuǎn)動(dòng)的方法和步驟如下：

(4)

式中：M為框架浮體的總力矩；F為第i單元上所受的合力；Ri為第i單元相對于質(zhì)心轉(zhuǎn)動(dòng)的力臂。

(5)

式中：I為浮架的轉(zhuǎn)動(dòng)慣量，aframe2表示浮架繞質(zhì)心轉(zhuǎn)動(dòng)的角加速度，wframe表示浮架繞質(zhì)心轉(zhuǎn)動(dòng)的角速度，w0表示浮架的初始角速度或上一時(shí)刻角速度，θ表示Δt時(shí)間內(nèi)浮架轉(zhuǎn)動(dòng)的角度。

根據(jù)布萊恩角計(jì)算公式[10-11, 16]可知，三維情況下Δt時(shí)間內(nèi)框架浮體轉(zhuǎn)動(dòng)的角度與該時(shí)刻的布萊恩角相等。具體浮架的轉(zhuǎn)動(dòng)和平動(dòng)的合成，即動(dòng)坐標(biāo)系和整體坐標(biāo)系的轉(zhuǎn)化公式如下：

(6)

式中：Vi為浮架第i單元整體坐標(biāo)系下的速度，由于每個(gè)單元Vframe和wframe均相等，所以無需使用i進(jìn)行標(biāo)記；P*i為浮架第i單元整體坐標(biāo)系下的x、y、z坐標(biāo)，依次為Pxi、Pyi、Pzi，P*center為浮架質(zhì)心位置整體坐標(biāo)系下該時(shí)刻的坐標(biāo)；θx、θy、θz表示x、y、z軸的布萊恩角。

1.2 網(wǎng)衣受力和運(yùn)動(dòng)

假定網(wǎng)衣是由有限的無質(zhì)量彈簧連接的集中質(zhì)量點(diǎn)所構(gòu)成，通過計(jì)算集中質(zhì)量點(diǎn)在波浪和浮架、錨繩作用下的位移，得到網(wǎng)衣的受力和形狀[16-18]。網(wǎng)衣的集中質(zhì)量點(diǎn)設(shè)于每個(gè)網(wǎng)目目腳的兩端，每個(gè)集中質(zhì)量點(diǎn)包含網(wǎng)目的1個(gè)網(wǎng)結(jié)和2個(gè)目腳，如圖2(b)所示，其中任意集中點(diǎn)4根網(wǎng)線的拉力沿順時(shí)針方向依次命名為f1～f4。關(guān)于單獨(dú)網(wǎng)衣的受力和運(yùn)動(dòng)可參見文獻(xiàn)[17-18]，其中給出了詳細(xì)的網(wǎng)目群化技術(shù)、集中質(zhì)量點(diǎn)法、網(wǎng)目局部坐標(biāo)建立等介紹，此處不再贅述。

下面主要介紹網(wǎng)衣網(wǎng)目恰好完整地系縛在浮架周圍上的處理方法以及系縛點(diǎn)位置網(wǎng)線拉力的計(jì)算方法。如圖2(b)為浮架沿圓形標(biāo)號1圓管展開后的結(jié)果，為了建模方便和使用一片網(wǎng)衣的思路，網(wǎng)衣的儲存變量形式為Pnet(ibar,1+4(jbar-1),3)，由于方形標(biāo)號和三角形標(biāo)號圓管連接端點(diǎn)有三處重復(fù)，所以網(wǎng)衣的列數(shù)應(yīng)為1+4(jbar-1)，ibar、4jbar的意義與上述一致。為了統(tǒng)一所有網(wǎng)目結(jié)點(diǎn)(質(zhì)量集中點(diǎn))均有相同數(shù)量拉力的網(wǎng)線，網(wǎng)衣系縛點(diǎn)并不集中質(zhì)量，而是集中在如圖2(b)實(shí)心圓點(diǎn)處，因此這些點(diǎn)中四角的端點(diǎn)集中了3.5倍的單根網(wǎng)線質(zhì)量以及其他點(diǎn)集中了3倍的單根網(wǎng)線質(zhì)量。與之對應(yīng)，平行于x和y軸的圓管上的系縛點(diǎn)受到相鄰網(wǎng)衣結(jié)點(diǎn)f1和f2的拉力，而平行于z軸的圓管上的系縛點(diǎn)受到相鄰網(wǎng)衣結(jié)點(diǎn)f1和f3、f2和f4的拉力。依次根據(jù)系縛點(diǎn)(浮架劃分單元點(diǎn))的坐標(biāo)和網(wǎng)目結(jié)點(diǎn)的坐標(biāo)分解網(wǎng)線拉力即可求出系縛點(diǎn)受力的矢量和。

如圖2(b)，平行于x和y軸的圓管的單元點(diǎn)編號和數(shù)量(jbar)與網(wǎng)衣不同，式(7)給出了它們之間的轉(zhuǎn)換關(guān)系，此關(guān)系是計(jì)算網(wǎng)衣和浮架的連接點(diǎn)(系縛點(diǎn))各物理量的關(guān)鍵。

(7)

式中：jp2和jp3依次對應(yīng)方形標(biāo)號和三角形標(biāo)號圓管的編號下標(biāo)up和down依次表示圖2(b)中上層圓管和下層圓管。j為任意網(wǎng)衣的列編號，其值為1～1+4(jbar-1)。

1.3 錨繩受力和運(yùn)動(dòng)

錨繩的數(shù)值模擬方法相對較成熟，采用集中質(zhì)量點(diǎn)法把錨繩劃分成假設(shè)用無質(zhì)量彈簧連接的球點(diǎn)可較好的描述錨繩的受力和運(yùn)動(dòng)，與彈簧受力不同的是錨繩僅在大于原長時(shí)才受力，如式(8)所示。錨繩集中質(zhì)量點(diǎn)受力和計(jì)算與網(wǎng)衣網(wǎng)目相似，主要包括重力、浮力、錨繩拉力、波浪拖曳力及波浪慣性力等[16-18]，具體計(jì)算公式見式(8)。

(8)

式中：Manchor表示錨繩質(zhì)量點(diǎn)的質(zhì)量，ΔManchor表示錨繩質(zhì)量點(diǎn)增加的慣性質(zhì)量；aanchor表示錨繩質(zhì)量點(diǎn)的運(yùn)動(dòng)加速度；Tanchor表示錨繩質(zhì)量點(diǎn)的三個(gè)方向的矢量拉力，而Tanchor表示相鄰錨繩質(zhì)量點(diǎn)因拉伸產(chǎn)生的標(biāo)量拉力值；FD-anchor、FI-anchor分別表示錨繩質(zhì)量點(diǎn)的拖曳力和慣性力；Ganchor表示錨繩質(zhì)量點(diǎn)的重力；Ffloat-anchor表示錨繩質(zhì)量點(diǎn)的浮力；l0為錨繩原始長度，l為變形后的長度；C1、C2為錨繩材料彈性系數(shù)；d為錨繩的直徑。

錨繩與浮架、錨繩與錨繩連接點(diǎn)的模擬方法需要特別處理，針對任一根錨繩連接點(diǎn)一定在錨繩首、尾質(zhì)量點(diǎn)處，模型建立時(shí)儲存數(shù)組仍保存錨繩首、尾質(zhì)量點(diǎn)的物理量，但某時(shí)刻內(nèi)僅計(jì)算一次錨繩首、尾質(zhì)量點(diǎn)，待該時(shí)刻所有其他質(zhì)量點(diǎn)和單元計(jì)算完成后，查找相應(yīng)連接點(diǎn)代替保存未計(jì)算錨繩首、尾質(zhì)量點(diǎn)的物理量(位置坐標(biāo)、速度等)。文中涉及到兩種浮架的錨泊方式，以圖3為例，一種系縛浮架的四角端點(diǎn)直接連接至海底，另一種先采用小直徑的輔助錨繩系縛浮架的四角端點(diǎn)，再連接一根大直徑的錨繩后至海底。第一種錨泊方式中，每根錨繩的首、尾點(diǎn)形式相同，尾點(diǎn)本身就固定于海底，首點(diǎn)與浮架的四角端點(diǎn)運(yùn)動(dòng)物理量重復(fù)，因此在模型中均不參與計(jì)算。待某時(shí)刻浮魚礁運(yùn)動(dòng)物理量全部計(jì)算完成后，把浮架的四角端點(diǎn)運(yùn)動(dòng)物理量直接賦值于錨繩的首點(diǎn)。第二種錨泊方式采用大、小兩種直徑的錨繩，4根輔助錨繩的尾點(diǎn)也是主錨繩的首點(diǎn)，因此在模型中4根輔助錨繩的首、尾點(diǎn)和主錨繩的尾點(diǎn)均不參與計(jì)算，單獨(dú)篩選出主錨繩的首點(diǎn)進(jìn)行計(jì)算，具體受力如圖3(b)所示。待某時(shí)刻浮魚礁運(yùn)動(dòng)物理量全部計(jì)算完成后，把浮架的四角端點(diǎn)運(yùn)動(dòng)物理量直接賦值于四根輔助錨繩首點(diǎn)，而把主錨繩的首點(diǎn)物理量賦值于輔助錨繩的尾點(diǎn)。這樣處理的好處主要體現(xiàn)在對每類浮魚礁構(gòu)成元素的描述均相同，能夠較方便地實(shí)現(xiàn)模型橫向拓?fù)浼昂筇幚怼?/p>

圖3 兩種錨泊方式計(jì)算方法Fig. 3 Calculation method of the two patterns of mooring

1.4 出水條件的判別方法

波浪作用下中上層浮魚礁的運(yùn)動(dòng)位置在波谷處時(shí)，浮魚礁各計(jì)算質(zhì)量點(diǎn)或單元可能發(fā)生出水情況。出水條件下將導(dǎo)致浮魚礁組成元素質(zhì)量點(diǎn)或單元不受浮力和波浪力，因此數(shù)值模型的建立需要進(jìn)行出水條件的判斷。判斷方法即浮魚礁各計(jì)算質(zhì)量點(diǎn)或單元的z坐標(biāo)與波面的關(guān)系，如式(9)所示。當(dāng)z大于η，水質(zhì)點(diǎn)速度和加速度均為0，相應(yīng)計(jì)算質(zhì)量點(diǎn)或單元所受浮力為0；當(dāng)z小于或等于η，考慮波浪水面變化的相對水深修正，如式(10)所示把z修正為zd，根據(jù)線性波理論[16]進(jìn)行計(jì)算。計(jì)算修正后的zd位置的波浪水質(zhì)點(diǎn)參數(shù)更符合出水條件下實(shí)際情況[18]。

=acos(ky-ωt)

(9)

(10)

式中：以靜水面位置為x軸，垂直向上為z軸，規(guī)定波浪沿y正方向傳播；η為水面變化，a為波浪振幅，k為波數(shù)，ω為圓頻率；dwave為水深，zd為修正后的位置坐標(biāo)。

特別指出，雖然出水條件下浮魚礁各計(jì)算質(zhì)量點(diǎn)或單元位置處的水質(zhì)點(diǎn)速度和水質(zhì)點(diǎn)加速度為0，但計(jì)算運(yùn)動(dòng)質(zhì)量點(diǎn)或單元拖曳力的莫里森公式中采用的是質(zhì)量點(diǎn)或單元相對速度(水質(zhì)點(diǎn)速度減質(zhì)量點(diǎn)或單元速度)，而該時(shí)刻的質(zhì)量點(diǎn)或單元速度一般不為0，需要進(jìn)一步地判斷，強(qiáng)制賦值為0。

2 數(shù)值模型驗(yàn)證

獲取驗(yàn)證資料的途徑主要有現(xiàn)場資料的采集和物理模型試驗(yàn)，由于文中設(shè)計(jì)的浮架式中上層浮魚礁尚未在海區(qū)現(xiàn)場投放，所以采用物理模型試驗(yàn)的方法驗(yàn)證數(shù)值模型的正確性。由于浮魚礁的運(yùn)動(dòng)狀態(tài)是其在波浪作用下各種力共同作用的結(jié)果，為了方便和簡化試驗(yàn)，選擇4根相同錨繩的錨泊方式，并且保證浮魚礁不出水，僅以浮魚礁上一點(diǎn)的運(yùn)動(dòng)變化驗(yàn)證數(shù)值模型的正確性。

2.1 物理模型制作

針對浮架的物理模型需要滿足幾何相似、動(dòng)力相似、運(yùn)動(dòng)相似，忽略圓管的剛度相似問題，即在原型中認(rèn)為圓管為剛性結(jié)構(gòu)。為保證浮架模型和原型之間重力和慣性力的相似，應(yīng)滿足弗勞德數(shù)Fr相等[19]。因?yàn)椴ɡ俗饔孟赂〖艿氖芰瓦\(yùn)動(dòng)帶有周期性變化，還應(yīng)保證周期性相似，即滿足斯特勞哈爾數(shù)相等，如式(11)所示。為保證網(wǎng)衣水動(dòng)力相似，物理模型制作中應(yīng)選擇網(wǎng)目較大的網(wǎng)衣，避免因網(wǎng)目過小導(dǎo)致的水流不相似問題。錨繩主要滿足彈性相似，可人工添加橡皮筋增加錨繩模型的彈性。

(11)

式中：V*、L*、T*為浮體結(jié)構(gòu)的特征速度、特性線尺度、周期，下標(biāo)m和p表示模型和原型。依據(jù)上述相似準(zhǔn)則，可推得原型與模型各物理量之間的比例關(guān)系，參見表1。

浮架使用PVC管、接頭和PVC膠黏接而成，除上下兩面以外，其他面周圍纏繞網(wǎng)衣，錨繩分別系縛在浮架的下底面端點(diǎn)。綜合考慮波浪水槽尺度32 m×0.8 m×1 m(長×寬×高)和實(shí)際海洋牧場海區(qū)面積，選取λ=20，建立物模模型的具體結(jié)構(gòu)參數(shù)見表2，實(shí)物如圖4。其中，表2中CD、CI、C1、C2等為計(jì)算浮架、網(wǎng)衣、錨繩波浪力的主要參數(shù)，具體取值參考Zhao等[9-10]，桂福坤[16]，陳天華等[17-18]的研究。

表1 原型與模型各參數(shù)間的比例關(guān)系Tab. 1 Proportional relationship of parameters between prototype and model

表2 物理模型制作參數(shù)Tab. 2 Physical model setting parameters

2.2 物理模型試驗(yàn)

物理模型試驗(yàn)在波浪水槽中進(jìn)行，水深0.75 m，波浪周期為1.1 s、1.4 s、2.1 s、4.0 s，對應(yīng)波浪高度為4.58 cm、16.63 cm、11.23 cm、18.36 cm，共四組波浪條件，如表3所示。波浪條件以圖4中浮魚礁模型上方的波高儀實(shí)測統(tǒng)計(jì)得出。利用CCD(charge coupled device，電荷耦合元件)相機(jī)采集圖4中示蹤點(diǎn)(發(fā)光二極管)連續(xù)的圖像，編寫程序?qū)D像進(jìn)行分析，獲得模型運(yùn)動(dòng)變化。圖像采集幀率會影響數(shù)據(jù)后期分析的準(zhǔn)確性，幀數(shù)越多，示蹤點(diǎn)運(yùn)動(dòng)軌跡越準(zhǔn)確，試驗(yàn)結(jié)果精度越高，幀數(shù)取10幀每秒即可滿足要求。因此，試驗(yàn)記錄時(shí)間為120 s，幀數(shù)為10幀每秒，每種波況平行重復(fù)做3次。

圖4 浮魚礁物理模型試驗(yàn)布置Fig. 4 Experiment layout of a floating fish attractor in a flume

開始造波前，浮魚礁處于靜止?fàn)顟B(tài)，選擇浮魚礁模型此時(shí)的位置作為參考標(biāo)準(zhǔn)，利用Matlab將不同時(shí)刻示蹤點(diǎn)的位置標(biāo)記在同一張圖片上，即得到示蹤點(diǎn)運(yùn)動(dòng)軌跡。使用游標(biāo)卡尺測量浮魚礁箱體管道之間的距離，同時(shí)計(jì)算浮魚礁在水中靜止?fàn)顟B(tài)該距離之間像素點(diǎn)個(gè)數(shù)，得到一個(gè)像素點(diǎn)空間中實(shí)際長度(距離/像素個(gè)數(shù))，在箱體不同位置重復(fù)多次計(jì)算，取平均值后一個(gè)像素在空間中表示長度為0.448 mm。先對CCD黑白圖像進(jìn)行灰度處理，得到各個(gè)像素點(diǎn)灰度值范圍在0～255之間的灰度圖片，其中黑色為0，白色為255。在白色點(diǎn)(255)周圍還有像素為254、253的像素點(diǎn)，文中只對像素灰度值為255的像素進(jìn)行計(jì)算[20]。

2.3 結(jié)果對比驗(yàn)證

以λ=20，按照表1中各參數(shù)間的比例關(guān)系建立浮魚礁數(shù)值模型，計(jì)算出示蹤點(diǎn)位置的運(yùn)動(dòng)軌跡后，再以λ=20縮小至物理模型尺度，物理模型試驗(yàn)和數(shù)值模擬的示蹤點(diǎn)位置軌跡對比見圖5，最大正、負(fù)向偏移和最大垂向偏移的對比見表3。選擇水深計(jì)算相對誤差率的原因是浮魚礁礁體的運(yùn)動(dòng)與錨繩的彈性形變相關(guān)。相同力條件下，錨繩越長，形變量越大，而水深決定著錨繩的長度。

表3 波浪要素及試驗(yàn)結(jié)果對比Tab. 3 Parameters of wave cases and experiment results analysis

圖5 不同工況下物理試驗(yàn)與數(shù)值模擬示蹤點(diǎn)運(yùn)動(dòng)軌跡對比Fig. 5 Movement tracks of the trace points of physical experiment and numerical modeling under four group waves

根據(jù)表3和圖5，對比物理模型試驗(yàn)和數(shù)值模擬結(jié)果發(fā)現(xiàn)，物理模型試驗(yàn)運(yùn)動(dòng)結(jié)果在水平方向上大于數(shù)值模擬結(jié)果；垂直方向上，左側(cè)部分驗(yàn)證較好，右側(cè)部分?jǐn)?shù)值模擬結(jié)果大于物理模型試驗(yàn)。以波浪條件的水深作為特征長度得到的水平方向和垂直方向的最大偏移相對誤差率在0.03～0.16之間。造成誤差的主要原因是水槽長時(shí)間造波存在反射現(xiàn)象，以及物理模型制作時(shí)還使用了細(xì)鐵絲系縛網(wǎng)衣，數(shù)值模型中并未考慮。綜合考慮，建立的數(shù)值模型能描述浮魚礁在波浪作用下的運(yùn)動(dòng)情況，也間接的驗(yàn)證了網(wǎng)衣和錨繩的受力。

3 中上層浮魚礁動(dòng)力響應(yīng)特性

中上層浮魚礁在波浪作用下可能為出水，因此根據(jù)第一節(jié)描述的數(shù)值方法計(jì)算了實(shí)際尺寸的浮魚礁運(yùn)動(dòng)和受力變化，其中錨泊為圖3中輔助錨和主錨相結(jié)合的方式。根據(jù)表2參數(shù)，浮魚礁數(shù)值模型尺寸以λ=20還原實(shí)際大小。浮架尺寸為1.69 m×1.69 m×4.8 m，圓管的直徑為0.2 m，密度相同；網(wǎng)衣的直徑0.003 m，目腳0.32 m；主錨繩直徑為0.02 m，輔助錨直徑為0.01 m，圖3中錨繩中間節(jié)點(diǎn)至浮架底部中心的距離為1 m；其余物理參數(shù)和計(jì)算參數(shù)均與表2中相同。根據(jù)舟山海域統(tǒng)計(jì)的5%累積頻率的波浪條件[21]為波高3 m，取水深12 m，波陡1/20，根據(jù)彌散方程計(jì)算的波周期為6.7 s，與文獻(xiàn)[21]統(tǒng)計(jì)的5%累積頻率的波周期6 s接近。

圖6給出了數(shù)值模型計(jì)算50 s后(浮魚礁運(yùn)動(dòng)已經(jīng)穩(wěn)定)的20 s內(nèi)浮架左右擺動(dòng)角度(波浪方向?yàn)檎⒉ɡ朔捶较驗(yàn)樨?fù))、錨繩拉力、網(wǎng)衣系縛點(diǎn)總拉力的歷時(shí)曲線(每0.1 s提取一個(gè)數(shù)據(jù)點(diǎn))?？梢缘贸?，三者存在明顯的周期性變化，最大錨繩拉力和最大系縛點(diǎn)總拉力均出現(xiàn)在擺動(dòng)角度突然變化的位置。如圖6中圓標(biāo)識的位置，擺動(dòng)角度負(fù)角度增大后驟然停止，致使錨繩拉力和系縛點(diǎn)總拉力達(dá)到最大且變化紊亂。圖6中矩形時(shí)刻為浮魚礁最大擺動(dòng)角度，此時(shí)錨繩拉力最小，且與系縛點(diǎn)總拉力變化趨勢相反。圖7為最大擺動(dòng)角度時(shí)刻的浮魚礁運(yùn)動(dòng)狀態(tài)及網(wǎng)衣受力分布，網(wǎng)衣受力僅在迎浪面和背浪面，并且拉力值遠(yuǎn)小于錨繩拉力。因此浮魚礁擺動(dòng)角度歷時(shí)曲線的光滑程度可以作為優(yōu)化浮魚礁結(jié)構(gòu)形狀和搭配浮體、配重的依據(jù)。

圖6 浮魚礁擺動(dòng)角度、錨繩拉力、網(wǎng)衣系縛點(diǎn)總拉力歷時(shí)曲線 Fig. 6 Time-varying of the floating fish attractor pendulum angle, the anchor rope tension and total tension of the netting fixing points

圖7 最大擺動(dòng)角度時(shí)刻的浮魚礁運(yùn)動(dòng)狀態(tài)及網(wǎng)衣受力分布Fig. 7 Kinestate of the floating fish attractor and netting force distribution at the moment of maximum pendulum angle

4 結(jié) 語

根據(jù)有限單元法和集中質(zhì)量點(diǎn)法，詳細(xì)闡述和推導(dǎo)了中上層浮魚礁各構(gòu)成組件受力和運(yùn)動(dòng)模擬計(jì)算方法，整體浮魚礁各組件連接位置處點(diǎn)或單元的相互耦合受力和運(yùn)動(dòng)處理方法，以及浮魚礁出水條件的判別方法和水質(zhì)點(diǎn)速度、加速度的修正方法。根據(jù)動(dòng)力相似和運(yùn)動(dòng)相似制作了浮魚礁的物理模型，并進(jìn)行了4組波浪條件的試驗(yàn)，對比分析浮架上一點(diǎn)的運(yùn)動(dòng)軌跡，驗(yàn)證了浮魚礁運(yùn)動(dòng)數(shù)值模型的正確性。最后建立了輔助錨和主錨相結(jié)合的錨泊方式的浮魚礁動(dòng)力響應(yīng)數(shù)值模型，通過分析浮魚礁擺動(dòng)角度、錨繩拉力、網(wǎng)衣系縛點(diǎn)總拉力歷時(shí)變化。主要結(jié)論總結(jié)如下：

1) 浮魚礁水動(dòng)力數(shù)值模型建立的關(guān)鍵是對其各組成要素的單元?jiǎng)澐趾涂臻g拓?fù)洌闹胁捎糜缮隙碌捻樞蚪o出了一種集浮架-網(wǎng)衣-錨繩為一體的拓?fù)浞绞?，以及浮魚礁運(yùn)動(dòng)和受力計(jì)算方法。該模型建立方法適用于小尺寸剛性和柔性結(jié)構(gòu)物及其組合而成的復(fù)雜結(jié)構(gòu)物的水動(dòng)力特性分析。

2) 單錨系泊條件下浮魚礁的錨繩拉力與上層礁體運(yùn)動(dòng)相關(guān)，擺動(dòng)角度驟然變化將導(dǎo)致錨繩拉力和系縛點(diǎn)總拉力達(dá)到最大值且變化紊亂。因此，浮魚礁擺動(dòng)角度歷時(shí)曲線的光滑程度可以作為優(yōu)化浮魚礁結(jié)構(gòu)形狀和搭配浮體、配重的依據(jù)。

3) 單錨系泊條件下浮魚礁網(wǎng)衣的受力面主要為迎浪面和背浪面，并且系縛點(diǎn)總拉力遠(yuǎn)小于錨繩拉力。因此，浮魚礁纏繞網(wǎng)衣在形成閉合空方體的同時(shí)，僅附加少量荷載，是浮魚礁設(shè)計(jì)的優(yōu)選形式。