三維纜索形態計算平臺的開發與應用

,,

(華中科技大學 船舶與海洋工程學院,武漢 430074)

1 纜索三維流場理論模型

基于以下假設：纜索不可深伸長；纜索只能承受拉力；不能承受彎矩及壓力；作用在纜索上的流體動力可以分解成切線方向的分力和法線方向的分力。設u，v，w為纜索上任意點P上的直角坐標單位矢量;u為P點處纜索方向。同樣設i，j和k為大地坐標系；dsu為P點處的長為無窮小的纜索元。φ為水平面與矢量u之間的夾角；φ為u在水平面的投影與矢量i的夾角，則

dx=dsu·i=dscosφcosφ

(1)

dy=dsu·j=dscosφsinφ

(2)

dz=dsu·k=dssinφ

(3)

設三維速度場(纜索相對于海水的速度)為V=Vxi+Vyj+Vzk，阻力的切向分力F和法向分力G、H的表達式就為

F=Fdsu=ρCdtπd(V·u)2dsu/2

(4)

G=Gdsv=ρCdnd(V·u)2dsv/2

(5)

H=Hdsw=ρCdnd(V·w)2dsw/2

(6)

式中：

V·u=Vxcosφcosφ+Vycosφsinφ+Vzsinφ

(7)

V·v=-Vxsinφ+Vycosφ

(8)

V·w=-Vxsinφcosφ-Vysinφsinφ+Vzcosφ

(9)

Cdt——纜索的切向阻力系數,Cdt=γCdn;

其中：Cdn——纜索的法向阻力系數。

纜索元上的作用力之和為零，平衡條件為

-Tu+(T+dT)(u+du)+Fdsu+Gdsv+
Hdsw-kpdsk=0

(10)

忽略高階無窮小量，得

(dT+Fds)u+(Tcosφdφ+Gds)v+
(Tdφ+Hds)w-pdsk=0

(11)

在u方向上：

dT+Fds-psinφds=0

(12)

在v方向上：

Tcosφdφ+Gds=0

(13)

在w方向上：

Tdφ+Hds-Pdscosφ=0

(14)

綜上可得纜索的平衡微分方程：

(15)

2 求解方法

上述微分方程組顯然是一階常微分方程組，根據特定的邊界條件對這些方程積分，得到纜索的平衡方程和拉力沿纜索徑跡的穩態分布。

本文采用變步長四階五級龍格-庫塔-芬爾格(Runge-Kutta-Fehlberg)方法[1-3]來求解。針對拖曳纜索的初值問題，即已知初始條件為

X0=[T0，φ0，φ0，x0,y0,z0]

將X0作為輸入條件沿著S=[0,S0]以一定步長進行迭代計算，其中S0為已知纜索總長度或者已知工作深度，就可以得出終端參數為

X1=[T1，φ1，φ1，x1,y1,z1]

同時步長點上的Xi也可以得出，這樣就可以得出整個拖曳纜索的形態。

3 GUI計算平臺

基于上述理論模型，利用Matlab圖形用戶界面功能，簡稱為GUI(graphic user interface)，編制計算平臺。

程序核心部分采用變步長四階五級龍格-庫塔-芬爾格方法，利用ODE模塊編寫計算模塊，同時利用GUI的I/O交互功能，對于纜索的一般初值問題，輸入初始參數即可求解。

平臺有文本框輸入初始參數，并顯示輸出結果，根據計算結果在圖像窗口自動繪制纜索三維形態。同時，實際工程中主要是已知繩長求深度和已知深度求繩長兩種問題，故平臺設置了工況選擇菜單，可以根據實際需要選擇不同工況求解。

4 實例

某拖曳器參數見表1。

表1 實例參數

4.1 無海流情況

如圖1，調用GUI界面進行計算，選擇“已知繩長求深度”工況，輸入初始參數，計算結果：纜索放出800 m時，潛水器的深度為304.7 m，此結果與原文獻結果十分接近。

4.2 有海流情況

有海流的情況，實質上是纜索相對于海水的速度矢量V=Vxi+Vyj+Vzk發生了變化，故求解方法與無海流情況無異。

對于有海流的情況，定義θ為海流速度矢量與拖曳速度矢量的夾角。當海流速度大小一定時，定性分析，θ從0～180°的變化過程中，海流阻礙拖曳的速度分量逐漸變大，那么必然使得拖曳纜索的深度逐漸變小。

圖1 GUI計算界面

假設存在與拖曳方向不一致的海流，其流速為0.5 kn；并且假定拖曳器受的阻力和凈浮力不變。當夾角為0、30、60、90、120、150及180°時，調用GUI界面進行計算，可以得出各角度下的纜索空間形態空間分布，見圖2，計算結果見表2。

圖2 不同海流夾角下纜索空間分布

表2 纜索深度隨海流夾角變化值

可以看出，在同一纜長和假定拖曳器受的阻力和凈浮力不變時，海流的存在對模型計算的影響是存在的，當海流速度為定值時，隨著夾角θ增大，纜索深度變小，也即意味著，在同一深度時，隨著夾角增大，所需纜長將增加。這與上述定性分析相吻合。

同時，可以對比θ與(180°-θ)的情況，事實上，這兩種情況，海流速度在拖曳方向的法向分量大小相同，故纜索的空間形態應在同一空間平面(曲面)內，只是由于在拖曳方向的分量方向相反，造成纜索的深度不同。

在不同角度下，纜索終端的拉力T也不同，故工程實際中可以根據計算結果選擇合適的拖曳方向，以在不影響拖曳目標的情況下實現拉力最小，從而達到節能目的。

4.3 考慮軸向變形的情況

上述理論模型假設纜索在拉力的作用下不存在軸向變形，下面進一步探討纜索在拉力的作用下伸長，存在軸向變形的情況

ds′=ds(1+ε)

(16)

(17)

(18)

式中：ε——應變；

ds′——拉伸后的微元長；

ds——拉伸前的微元長；

A——拉伸后纜索的橫截面；

d0——拉伸前纜索微元段的直徑；

d1——拉伸后纜索微元段的直徑；

E——材料的彈性模量。

根據式(16)～(18)，可推導出應變ε的表達式為

(19)

于是考慮彈性后的纜索三維模型微分方程組為

(20)

根據上述方程組，對上例平臺中的方程組做出相應修改，調用GUI界面計算，對于銅制(或黃銅，青銅)纜芯，E=70～130 GPa，此例取為70 GPa，計算結果見圖3。

圖3 考慮彈性影響下的纜索曲線

可以看出，考慮彈性和不考慮彈性的情況下，兩條曲線的吻合程度非常好，這是由于對于纜芯為銅制的纜索，彈性模量E一般較大，故應變就很小，所以對于常用纜索材料，彈性對于纜索形態的影響可以忽略不計。

[1] RECKTENWALD G.數值分析和MATLAB實現與應用[M].伍衛國,萬 群,張 輝,等.譯.北京：機械工業出版社，2004.

[2] 薛定宇.高等應用數學問題的Matlab求解[M].北京：清華大學出版社,2008.

[3] 張圣君.浮標拖曳纜索形狀及其力學計算[D].武漢：華中科技大學，2003.