水下多纜多體拖曳系統(tǒng)運(yùn)動(dòng)建模與模擬計(jì)算

王 飛，丁 偉，鄧德衡，吳小峰

(1. 上海交通大學(xué) 海洋工程國家重點(diǎn)實(shí)驗(yàn)室； 高新船舶與深海開發(fā)裝備協(xié)同創(chuàng)新中心， 上海 200240； 2. 武漢第二船舶設(shè)計(jì)研究所， 武漢 430064； 3. 中國船舶科學(xué)研究中心， 江蘇 無錫 214082)

水下拖曳系統(tǒng)作為一種基本的水下探測裝備，廣泛應(yīng)用于海底探測、軍事反潛、海工應(yīng)用等領(lǐng)域，能否準(zhǔn)確預(yù)知其在各種情況下的空間姿態(tài)及動(dòng)力響應(yīng)特性，是拖曳系統(tǒng)設(shè)計(jì)與應(yīng)用中的重要問題.拖曳系統(tǒng)通常由水面拖曳母船、牽引纜、水下拖體、陣列、尾繩等設(shè)備組成，拖曳形式也因用途的不同而千變?nèi)f化[1-2].現(xiàn)代新型拖曳系統(tǒng)為了增強(qiáng)探測能力，開始向多陣列(拖纜)多拖體系統(tǒng)發(fā)展，一個(gè)系統(tǒng)中可包含數(shù)個(gè)陣列、數(shù)個(gè)拖體或其組合，效率顯著提升；但系統(tǒng)也因此變得更為復(fù)雜，實(shí)際應(yīng)用中也暴露出一些問題，如運(yùn)動(dòng)響應(yīng)特性不明晰、纜索糾纏、海底觸碰，這些問題常會對系統(tǒng)的正常運(yùn)行造成影響.

水下拖曳系統(tǒng)動(dòng)力學(xué)的理論研究在很早就已經(jīng)廣泛展開，國內(nèi)外眾多學(xué)者針對各自的研究目標(biāo)，建立運(yùn)動(dòng)數(shù)學(xué)模型，采用數(shù)值模擬等方法來研究其在不同情況下的運(yùn)動(dòng)響應(yīng)特性，包括各種穩(wěn)態(tài)運(yùn)動(dòng)研究及動(dòng)態(tài)運(yùn)動(dòng)研究.穩(wěn)態(tài)運(yùn)動(dòng)研究方面，Wang等[3]在前人研究結(jié)果的基礎(chǔ)上，針對多種形式的拖曳系統(tǒng)，較為全面地考慮了各種情況下的穩(wěn)態(tài)運(yùn)動(dòng)，建立初始值方程，提出了多維嵌套二分法進(jìn)行求解，以解決諸多實(shí)際問題，并將其應(yīng)用到系統(tǒng)設(shè)計(jì)中去.系統(tǒng)的動(dòng)態(tài)運(yùn)動(dòng)研究更為廣泛，通常研究人員在拖纜運(yùn)動(dòng)微分控制方程的基礎(chǔ)上采用集中質(zhì)量法、有限差分法等，來分析其在各種情況下的動(dòng)態(tài)運(yùn)動(dòng)[4-12].Huang[4]最早提出了集中質(zhì)量法來計(jì)算拖纜的動(dòng)態(tài)運(yùn)動(dòng)響應(yīng)，它們能夠很好地預(yù)報(bào)拖纜系統(tǒng)在正常狀態(tài)下的運(yùn)動(dòng).在此基礎(chǔ)上，國內(nèi)外眾多學(xué)者展開了卓有成效的拓展性研究[5-12]，內(nèi)容涉及水下拖體耦合動(dòng)力學(xué)研究[5,7-9]、拖曳系統(tǒng)導(dǎo)桿影響分析[10]、系統(tǒng)收放模擬[11]、海底纜鋪設(shè)運(yùn)動(dòng)研究[12].

多陣列拖曳系統(tǒng)作為近些年出現(xiàn)的一種系統(tǒng)，相關(guān)研究文獻(xiàn)較少，其在不同情況下的運(yùn)動(dòng)響應(yīng)特性不甚明晰.為此本文針對水下多陣列多拖體拖曳系統(tǒng)，展開基礎(chǔ)研究工作，探討其在不同情況下的穩(wěn)態(tài)和動(dòng)態(tài)運(yùn)動(dòng)響應(yīng)特性.擬采用集中質(zhì)量法建模思想，結(jié)合潛艇6自由度運(yùn)動(dòng)控制方程，建立拖纜(陣列)及拖體的運(yùn)動(dòng)控制方程，并通過建立不同的連接耦合邊界條件將其耦合成一個(gè)整體.通過數(shù)值模擬研究系統(tǒng)在各種不同情況下的運(yùn)動(dòng)響應(yīng)特性，以期給出一些規(guī)律性結(jié)論.

1 系統(tǒng)運(yùn)動(dòng)建模

本文將拖纜(陣列)視為理想的柔性纜索，忽略纜彎曲及扭轉(zhuǎn)剛度的影響； 鑒于水面拖曳母船的排水量逾3 000 t，且配有航跡姿態(tài)控制系統(tǒng)，本文忽略纜對水面拖曳母船的耦合作用；因各個(gè)纜的運(yùn)動(dòng)控制方程均一致，簡便起見A、B、C、D 4根拖纜(陣列)的參數(shù)符號全文統(tǒng)一表示，不再單獨(dú)指明.

典型的多陣列拖曳系統(tǒng)如圖1所示，其中水下拖體內(nèi)部載有多種測控設(shè)備，并連接4根拖纜，纜A為牽引纜，水中長度可調(diào)整用以進(jìn)行變深控制，纜B和纜C為近乎0浮力的陣列纜，纜D下端為一小型拖曳設(shè)備.水下拖體尺度較大，對系統(tǒng)運(yùn)動(dòng)有重要影響，本文采用6自由度運(yùn)動(dòng)建模，而纜D下端的小型設(shè)備則將其視為質(zhì)量點(diǎn)進(jìn)行簡化建模處理.對于每一根拖纜(陣列)，s為纜長坐標(biāo)，S為其總長，其中拖纜下端/尾端的坐標(biāo)記為s=0，(θ,φ)為拖纜微元姿態(tài)角，方位角θ為纜順時(shí)針偏離y軸的角度，φ為抬升角，這兩個(gè)角度取決于纜的空間位置.水下拖體的3個(gè)姿態(tài)角，參考潛艇操縱性運(yùn)動(dòng)分別定義為：艏向角ψ，為艏向水平偏離y軸的角度，順時(shí)針為正；縱傾角?，抬艏為正；橫傾角φ，右傾為正.

圖1 水下多纜多體拖曳系統(tǒng)Fig.1 Multi-cable multi-body towed system

1.1 坐標(biāo)系統(tǒng)

建立3個(gè)坐標(biāo)系統(tǒng)描述整個(gè)系統(tǒng)的運(yùn)動(dòng)：慣性坐標(biāo)系(Oxyz)，拖纜局部坐標(biāo)系(Obtn)和拖體/拖船運(yùn)動(dòng)坐標(biāo)系(O′ξηζ).其中慣性坐標(biāo)系原點(diǎn)位于水面，z軸垂直向上；拖纜局部坐標(biāo)系附于拖纜上，t軸為切向，n為法向，b為側(cè)法向；拖體運(yùn)動(dòng)坐標(biāo)系隨拖體運(yùn)動(dòng)，ξ軸指向船首，η軸指向右舷，ζ軸指向下方.這3個(gè)坐標(biāo)系通過姿態(tài)角相互關(guān)聯(lián)，其轉(zhuǎn)換關(guān)系為

(1)

式中：Q與A分別為拖體和拖纜坐標(biāo)系轉(zhuǎn)換矩陣.

1.2 拖纜運(yùn)動(dòng)控制方程

拖纜微元的控制方程為

(2)

式中：M為拖纜微元質(zhì)量矩陣；T為張力；D為阻力；G為重力；B為浮力；x為其空間位置向量.

根據(jù)集中質(zhì)量法思想，將一拖纜從尾端至上端離散為N段，即N+1個(gè)節(jié)點(diǎn)，其中：尾端s=0，為第i=0個(gè)節(jié)點(diǎn)；上端s=S，為第i=N個(gè)節(jié)點(diǎn).

對第i個(gè)節(jié)點(diǎn)應(yīng)用牛頓第二定律，得拖纜節(jié)點(diǎn)的基本運(yùn)動(dòng)控制方程：

(3)

1.3 拖體運(yùn)動(dòng)控制方程

XH+XS+∑Xci

YH+YS+∑Yci

ZH+ZS+∑Zci

MH+MS+∑Mci

其中：(u,v,w,p,q,r)為拖體的線速度及角速度；I為各個(gè)方向上的轉(zhuǎn)動(dòng)慣量；方程右側(cè)(X,Y,Z,K,M,N)為作用在拖體上的力和力矩；下標(biāo)H表示拖體水動(dòng)力；下標(biāo)S表示拖體重浮力作用分量；下標(biāo)c表示拖纜作用力.拖纜作用力將作為邊界條件在1.4節(jié)給出.

拖體的重浮力作用力為

作為運(yùn)動(dòng)關(guān)聯(lián)條件，水下拖體的姿態(tài)角同角速度有如下關(guān)系：

1.4 邊界條件

在拖曳系統(tǒng)運(yùn)動(dòng)過程中，拖纜與水下拖體/水面拖船是相互作用、相互影響的.根據(jù)不同的情況，有3種不同的邊界條件：拖纜自由端邊界條件、拖纜下端邊界條件和水下拖體(水面拖船)-纜耦合條件.為簡便起見，一些公式中將拖纜節(jié)點(diǎn)的下標(biāo)忽略不寫.

(1) 拖纜自由端邊界條件.對于沒有拖體的自由端(纜B和纜C)，視其為一節(jié)點(diǎn)，應(yīng)用集中質(zhì)量法控制方程：

(4)

式中：I為單位矩陣；l為節(jié)點(diǎn)間長度；g為重力加速度；下標(biāo)1/2表示節(jié)點(diǎn)0和1之間的物理量.

(2) 拖纜下端小型設(shè)備邊界條件(纜D).該處設(shè)備尺度相對很小，可直接將其溶入到纜的控制方程，采用3自由度運(yùn)動(dòng)方程進(jìn)行求解：

(5)

式中：Md為小型設(shè)備的3自由度質(zhì)量及附加質(zhì)量矩陣；Fd為其上的作用力，

右式第1項(xiàng)為其重浮力，第2項(xiàng)為其運(yùn)動(dòng)阻力，Sd為設(shè)備的特征面積，Cd為其阻力系數(shù).

(3) 水下拖體(水面拖船)-纜耦合邊界條件.在拖纜與拖體連接點(diǎn)處邊界條件包括運(yùn)動(dòng)條件和動(dòng)力條件，本文將連接點(diǎn)處拖體的位置和速度作為運(yùn)動(dòng)耦合邊界條件，以確定拖纜連接節(jié)點(diǎn)的位置速度；而拖纜的張力作為動(dòng)力耦合邊界條件，融入到拖體操縱運(yùn)動(dòng)控制方程.

(6)

式中：t′為時(shí)間；(xt,yt,zt)表示拖體在慣性坐標(biāo)系下的位置；Rtp為拖體坐標(biāo)系下連接點(diǎn)的坐標(biāo)值.

作為動(dòng)力耦合邊界條件，拖纜作用于拖體的作用力(如纜A下端，其余類似)，經(jīng)過坐標(biāo)轉(zhuǎn)換可以表示為

水面拖船的位置與速度為整個(gè)系統(tǒng)的輸入，是已知的時(shí)間函數(shù)；而水下拖體的位置與速度則由拖體的運(yùn)動(dòng)控制方程計(jì)算得到.

1.5 系統(tǒng)穩(wěn)態(tài)運(yùn)動(dòng)求解

系統(tǒng)的穩(wěn)態(tài)解通常作為動(dòng)態(tài)運(yùn)動(dòng)的初始條件.各拖纜的穩(wěn)態(tài)運(yùn)動(dòng)控制方程在局部坐標(biāo)系下為

(7)

式中：u、d、σ分別為纜的單位長度質(zhì)量、直徑、截面積；u=(ub,ut,un)為局部坐標(biāo)系下纜的速度；T為拖纜張力；Ct、Cn分別為切向、法向阻力系數(shù)；ε為纜的應(yīng)變.

拖纜位置在慣性系下的微分關(guān)系為

對此穩(wěn)態(tài)運(yùn)動(dòng)控制方程，采用4階龍格庫塔積分從拖纜下端開始沿纜長積分，得到整個(gè)拖纜的穩(wěn)態(tài)解，拖纜在下端的初始值(T0,θ0,φ0)根據(jù)不同的情況分別如下：

(1) 拖纜自由端初始值.對于纜B和纜C的自由端，拖纜張力與歐拉角的變化率均為0.因拖纜側(cè)向無作用力，姿態(tài)角

θ≡ const=ψ

所以可以將穩(wěn)態(tài)問題轉(zhuǎn)換到二維空間下求解，最終拖纜的初始值為

(8)

式中：

ww=(μ-ρσ)g

ww為纜在水中單位長度重量，

(2) 水下拖體連接點(diǎn)初始值.對于纜A以及纜D的下端，可以由作用在拖體上的外力來確定拖纜的初始值，

(9)

式中：Fx、Fy、Fz為作用在水下拖體上除該拖纜以外的外力.

(3) 水下拖體穩(wěn)態(tài)求解.由于拖體的控制方程是一個(gè)高階的非線性系統(tǒng)，其穩(wěn)態(tài)解不能解析計(jì)算，本文在纜B以及纜C張力穩(wěn)態(tài)解的基礎(chǔ)上，采用動(dòng)態(tài)計(jì)算足夠長時(shí)間的方式得到其穩(wěn)態(tài)解(u,w,?)，而拖體的橫向速度和橫傾角均為0，艏向角則為系統(tǒng)的艏向角.

2 控制方程數(shù)值求解

由式(3)可得各拖纜的運(yùn)動(dòng)控制方程：

(10)

將式(10)與拖體運(yùn)動(dòng)控制方程和邊界條件聯(lián)立，可得到整個(gè)系統(tǒng)的運(yùn)動(dòng)方程，并采用4階龍格庫塔方法在時(shí)域內(nèi)積分求解.由于控制方程是條件穩(wěn)定的，計(jì)算時(shí)間步長必須很小，一般為毫秒級.

3 仿真計(jì)算與分析

整個(gè)拖曳系統(tǒng)由4個(gè)拖纜、1個(gè)拖體和1個(gè)小型拖曳設(shè)備組成，屬短陣高航速系統(tǒng)，其主要參數(shù)如表1和表2所示，海水密度取為ρ=1 025 kg/m3，無海流，牽引纜上端連接于水面拖曳母船，水面母船的位置與速度為整個(gè)系統(tǒng)的輸入量.4個(gè)拖纜同拖體連接點(diǎn)坐標(biāo)分別為A(0, 0, 0) m，B(-1.7, -1.1, 0) m，C(-1.7,1.1,0) m，D(-0.2,0,0.62) m.整體上拖纜以5 m的間距進(jìn)行離散，在曲率變化較大的地方適當(dāng)加密，拖纜節(jié)點(diǎn)數(shù)量共126個(gè)，計(jì)算時(shí)間步長取2 ms.

表1 各拖纜物理參數(shù)Tab.1 Physical parameters of the towed cables

表2 水下拖體參數(shù)及無因次慣性類水動(dòng)力導(dǎo)數(shù)

Tab.2 Parameters and non-dimensional inertial derivatives of the underwater towed body

參數(shù)取值長度/m2.72質(zhì)量/kg673體積/m30.312重心/m(0.02, 0, 0.46)浮心/m(0.00, 0, 0.22)Z'q·,M'w·0.01570K'p·-0.00207K'r·,N'p·0.00020M'q·-0.00278N'r·-0.00371X'u·-0.0067Y'v·-0.0763Y'p·,K'v·0.00682Y'r·,N'v·0.0159Z'w·-0.0861Ixx/(kg·m2)170.4Iyy/(kg·m2)197.6Izz/(kg·m2)175.3Ixz/(kg·m2)6.2

鑒于水下拖體一般無大幅機(jī)動(dòng)，水動(dòng)力僅取主要項(xiàng)，各無因次導(dǎo)數(shù)分別為

纜D下端小型拖曳設(shè)備的參數(shù)為：質(zhì)量61 kg，體積0.028 m3，其質(zhì)量矩陣(單位：kg)及作用力(單位：N)近似為

3.1 穩(wěn)態(tài)運(yùn)動(dòng)計(jì)算與分析

系統(tǒng)的穩(wěn)態(tài)計(jì)算結(jié)果作為下文動(dòng)態(tài)運(yùn)動(dòng)計(jì)算的初始條件.

本節(jié)拖曳速度v=4,6,8,10 kn，方向沿y軸正方向，水面拖船拖曳點(diǎn)位置為(0,0,0) m.經(jīng)計(jì)算得到不同速度時(shí)系統(tǒng)的穩(wěn)態(tài)結(jié)果， 如表3、圖2和圖3所示.其中牽引纜上端的張力分別為Ti=N-1/2=5.7,8.6,13.5,20.2 kN，陣列纜尾端在水下的深度分別為z=-74.2,-47.6,-33.7,-25.2 m.

表3 拖體的穩(wěn)態(tài)參數(shù)Tab.3 Steady state parameters of the towed body

圖2 穩(wěn)態(tài)纜形Fig.2 Steady state configurations

圖3 穩(wěn)態(tài)張力分布Fig.3 Steady state tension distribution

結(jié)果顯示，陣列纜尾端較前端稍有些上浮，因?yàn)殛嚵欣|的密度低于水密度，有正浮力.而水下拖體縱傾角由于受到拖體重浮力力矩、水動(dòng)力縱傾力矩和纜張力形成的力矩多種因素的影響，所以不同速度時(shí)變化較復(fù)雜，低速時(shí)為負(fù)值，高速時(shí)逐漸增大為正值，且不具有線性關(guān)系.

3.2 加速運(yùn)動(dòng)計(jì)算

加速運(yùn)動(dòng)由4 kn加速至10 kn，加速時(shí)長分別取為100和400 s.圖4所示為牽引纜上端張力變化曲線，圖5所示為系統(tǒng)各部位深度變化曲線.由圖可見，纜上張力不論加速時(shí)間長短，加速過程中逐步增加，在加速過程結(jié)束后很快達(dá)到穩(wěn)態(tài)，沒有出現(xiàn)超調(diào)；而拖體和陣列纜在水下的深度，則沒能很快地達(dá)到穩(wěn)態(tài)，加速時(shí)間越短達(dá)到穩(wěn)態(tài)的時(shí)間越長，特別是零浮力的陣列纜.

圖4 加速運(yùn)動(dòng)牽引纜上端張力變化曲線Fig.4 Tension at the lead-in cable upper end for acceleration motion

圖5 加速運(yùn)動(dòng)系統(tǒng)深度變化曲線Fig.5 Depth of the towed system for acceleration motion

3.3 回轉(zhuǎn)運(yùn)動(dòng)計(jì)算

拖曳系統(tǒng)作360° 勻速回轉(zhuǎn)，然后繼續(xù)沿原方向直線航行，計(jì)算參數(shù)為：回轉(zhuǎn)半徑取R=100,200 m；回轉(zhuǎn)速度取v=10 kn.圖6所示為回轉(zhuǎn)過程中系統(tǒng)的響應(yīng)曲線，包括拖體深度、拖體橫傾角、縱傾角及纜張力的變化， 圖7所示為各參考位置的軌跡歷時(shí)曲線.

結(jié)果顯示：兩個(gè)水聽器陣列纜盡管間距很小，但回轉(zhuǎn)中始終接近平行運(yùn)行，無交叉纏繞現(xiàn)象.拖體回轉(zhuǎn)過程中在水下的深度有較大幅度變化，特別是小半徑回轉(zhuǎn)時(shí)，深度由初始時(shí)26 m降低至47.6 m，增加了83%，如水深有限則可能發(fā)生碰撞事故；與此同時(shí)纜上的張力大幅降低，拖體的橫、縱傾角也有較大幅度變化，且有一定的沖擊變化現(xiàn)象； 同時(shí)系統(tǒng)上的各參考位置，其回轉(zhuǎn)半徑均小于拖船的回轉(zhuǎn)半徑，特別是陣列纜尾部，回轉(zhuǎn)半徑越小減小幅度越厲害，此時(shí)陣列纜有可能和牽引纜碰撞或纏繞，導(dǎo)致系統(tǒng)不能正常工作.經(jīng)仔細(xì)計(jì)算觀測，在R=100 m，t=144 s 左右時(shí)陣列纜末端同牽引纜之間的點(diǎn)線距離已接近于零，碰撞纏繞現(xiàn)象可能已發(fā)生.

圖6 拖曳系統(tǒng)回轉(zhuǎn)運(yùn)動(dòng)響應(yīng)曲線Fig.6 Dynamic response of the towed system during loop maneuver

圖7 系統(tǒng)回轉(zhuǎn)軌跡曲線Fig.7 Trajectory of the towed system during loop maneuver

3.4 水面拖船振蕩響應(yīng)計(jì)算

分別考慮3個(gè)方向上拖船振蕩(橫蕩、縱蕩、垂蕩)對于系統(tǒng)的影響，幅值取為0.5 m，周期為10 s；拖曳速度為v=4,10 kn.圖8～10給出了穩(wěn)定后(t′>800 s)，系統(tǒng)各參考位置處偏蕩歷時(shí)曲線，圖中Δx,Δy,Δz為各個(gè)方向上的偏蕩變化量； 其中橫蕩時(shí)拖船處(纜A上端)張力由于幾乎沒有變化而沒有在圖中給出，纜C的結(jié)果同纜B全都類似，也沒有給出.

圖8 橫蕩運(yùn)動(dòng)響應(yīng)Fig.8 Response for sway motion

圖9 縱蕩運(yùn)動(dòng)響應(yīng)Fig.9 Response for surge motion

圖10 垂蕩運(yùn)動(dòng)響應(yīng)Fig.10 Response for heave motion

結(jié)果表明，拖船的橫蕩運(yùn)動(dòng)對于系統(tǒng)的影響很小，拖纜張力幾乎沒有變化，各參考位置的橫蕩幅度也極小，可以忽略.但是縱蕩和垂蕩運(yùn)動(dòng)對于整個(gè)系統(tǒng)有著顯著的影響，拖纜張力及各參考位置的深度均有較大幅度的波動(dòng)，同時(shí)結(jié)果顯示，隨著拖曳速度的增加，縱蕩和垂蕩運(yùn)動(dòng)的影響程度會逐漸減小(參考位置偏蕩幅度減小)，但橫蕩的影響會稍有增大.

4 結(jié)語

本文采用數(shù)值仿真方法，針對水下多纜多體拖曳系統(tǒng)，通過建立其運(yùn)動(dòng)數(shù)學(xué)模型，展開了較為系統(tǒng)的研究，以探討其在不同情況下的運(yùn)動(dòng)響應(yīng)特性.其中水下拖纜(陣列)采用集中質(zhì)量法進(jìn)行建模，水下拖體采用潛艇6自由度運(yùn)動(dòng)控制方程，通過建立拖纜和拖體耦合條件以形成整個(gè)系統(tǒng)的運(yùn)動(dòng)數(shù)學(xué)模型，采用4階龍格庫塔方法積分求解.最后針對4纜2拖體拖曳系統(tǒng)展開數(shù)值模擬研究，分析其在不同情況下的穩(wěn)態(tài)和動(dòng)態(tài)運(yùn)動(dòng)響應(yīng)特性，給出了一些規(guī)律性結(jié)論：

(1) 系統(tǒng)回轉(zhuǎn)時(shí)，拖曳深度和纜上張力有大幅度地降低，拖體的橫、縱傾角在回轉(zhuǎn)結(jié)束時(shí)會出現(xiàn)較大的沖擊現(xiàn)象，同時(shí)發(fā)現(xiàn)高速小半徑回轉(zhuǎn)為最危險(xiǎn)，拖體及陣列纜旋回半徑大幅減小，陣列纜可能同牽引纜發(fā)生碰撞纏繞現(xiàn)象.

(2) 水面拖船橫蕩運(yùn)動(dòng)對于系統(tǒng)的影響很小，而縱蕩和垂蕩卻對拖曳深度及纜張力有顯著影響，其影響程度負(fù)相關(guān)于航速的變化.