左衛(wèi)廣，王永學

（大連理工大學 海岸和近海工程國家重點實驗室，大連116024）

沉管法是當前海底隧道建造施工的新方式，而沉管管段的安全沉放是沉管隧道施工的關鍵環(huán)節(jié)。駁船吊掛沉放法作為一種常用的沉放方式，其沉管管段的運動特性是波浪直接引起的沉管管段運動與波浪作用下駁船施工裝備的運動間接引起的沉管管段運動的復合運動。駁船升沉運動和波浪聯合作用下的沉管管段的沉放運動特性是一個十分復雜的水動力學問題。

目前國內外關于沉管管段的沉放運動特性研究工作大都是忽略了駁船運動的影響。詹德新［1］試驗研究了沉管沉放過程中的流體作用力和纜繩纜力，考慮了流速和波高等環(huán)境因素的影響，并分析了壓載水對沉管沉放穩(wěn)定性的影響。Jensen［2］等建立了施工海域的波浪模型，分析了沉管的沉放運動特性，研究了近海波況對沉管隧道施工沉放的影響。陳智杰［3］對沉管管段在固定深度下的運動特性進行了試驗研究，并采用邊界元方法對沉管管段沉放運動進行了時域數值模擬，數值模擬結果與試驗結果比較，吻合良好。

實際工程中，當沉管管段進行沉放施工時，沉管管段不僅受到波浪和水流的作用，還受到海上駁船運動的影響。駁船在波浪和水流作用下的運動，通過纜繩的傳遞將會對沉管運動有很大影響。本文主要研究駁船升沉運動和波浪聯合作用對沉管管段運動特性的影響。建立了駁船升沉運動和波浪聯合作用下沉管管段運動的數學模型，通過數值算例分析，探討了駁船做不同振幅與不同相位升沉運動情形下的沉管管段運動特性。

1 波浪作用下沉管運動數學模型

駁船吊掛沉放法進行沉管管段沉放作業(yè)的示意圖見圖1。在圖1 中，引入2 個坐標系統(tǒng)：參考坐標系O0-x0y0z0和動坐標系O-xyz。坐標平面O0-x0y0z0位于靜水面，x0軸為沿沉管的長度方向，z0軸為垂直水面向上；坐標系O-xyz 與管段固結，隨管段一起搖蕩，x 軸為沿管段的長度方向。當管段處于平衡位置時，動坐標系Oxyz 和參考坐標系O0-x0y0z0是重合的。參考坐標系是表征沉管管段搖蕩位移和姿態(tài)的基準，當管段由于外力的作用產生偏移時，偏移量由動坐標系相對于參考坐標系的運動來描述。

圖1 所示的沉管管段的運動可以看作是兩部分運動的合成，一部分是忽略駁船運動情形（假定吊纜上端點為固定情形）下波浪作用引起的沉管管段運動；另一部分是由駁船運動引起的沉管管段運動。

圖1 沉管管段示意圖Fig.1 Sketch of immersing tunnel element

本文采用直接時域法求解波浪作用引起的沉管管段運動?；趧萘骼碚?，流場的速度勢滿足Laplace 方程的定解問題。應用格林定理建立時域內的積分方程，使用邊界元方法（Hess 和Smith［4］；Garrison［5］；Liu 等［6］；Teng 和Eatock Taylor［7］）數值求解時域內的積分方程得到速度勢后，可計算出沉管管段所受到的波浪力，求解沉管管段的運動方程得到波浪作用引起的沉管管段運動響應。

1.1 流場控制方程

流場內，速度勢φ 滿足Laplace 方程

式中：A 為波幅；β 為浪向角；g 為重力加速度；h 為水深；k 為波數；ω 為入射波的頻率；t 為時間；φm為φ，φI或φS中的任何一個；p 為流場中點的位置矢量；n 為物體表面的法向量，物面向外為正；Vn為物體表面上的法向速度。

1.2 速度勢求解

本文采用格林函數法來處理速度勢的求解問題。有限水深條件下時域格林函數可以表達成以下兩式之和的形式［8－9］

G0和Gf分別稱為瞬時項和記憶項，其定義為

式中：τ 為時間；p（x，y，z）和q（ξ，η，ζ）分別為場點和源點矢量；J0為零階貝塞爾函數；r 為場點和源點的距離；R 為場點和源點的水平距離；r2為場點和源點關于水底面鏡像的距離。

對速度勢的初邊值問題應用格林定理，考慮線性問題，物面邊界用平均邊界SB代替SB（τ），經過推導可得到在自由面SF（τ）、物面SB、海底邊界SD和無窮遠邊界S∞組成的流體域Ω 內的時域內積分方程

應用邊界元方法，數值求解時域內積分方程（10），可得到流場速度勢φS。

1.3 荷載計算

沉管在沉放過程中，受到波浪力和纜繩纜力的共同作用。

波浪力可表示為

式中：ρ 為流體密度；速度勢φ 包括入射勢φI和散射勢φS。

纜繩纜力可表示為

式中：T 為纜繩張力；Ce 為纜繩彈性系數；d 為纜繩直徑；ΔS 為纜繩變形長度；S 為纜繩原始長度；n 為指數。

當纜繩為鋼纜時，n＝3/2。

1.4 沉管運動方程

沉管運動方程可寫為如下形式

式中：Fk為總的廣義水動力荷載分量；Gk為纜繩沉放系統(tǒng)對沉管施加的外部作用力和力矩；M 為質量矩陣；B為系統(tǒng)阻尼矩陣；C 為恢復力矩陣，由于沉管位于水面以下，故靜恢復力為0；ξj為第j 個模態(tài)下沉管管段的響應幅值，j=1～6 分別為沉管的縱蕩、橫蕩、升沉、橫搖、縱搖和回轉運動。

沉管運動方程式可表示為如下的二階微分方程形式

應用四階Runge-Kutta 方法求解式（14），可得到沉管管段的位移和速度表達為

式中：M1、M2、M3和M4分別為

2 考慮駁船升沉運動情形下沉管運動數學模型

假定駁船在升沉方向做強迫簡諧運動，其升沉運動方程采用如下表達式

式中：b 為駁船做升沉運動的振幅；ω 為入射波頻率；θ 為駁船運動與波浪之間的相位差；tm為緩沖周期。

在t 時刻，駁船的位置為z（t），同時考慮波浪作用與駁船升沉運動影響的沉管位置為ξj（t）（j=1，2，3，4，5，6）。在t+Δt 時刻，首先考慮駁船升沉運動對沉管位置的改變。駁船在t+Δt 時刻的位置z（t+Δt），相對于t 時刻的瞬時位置變化為

在僅考慮駁船升沉運動對沉管運動的影響時，認為駁船做升沉運動瞬時纜繩長度不發(fā)生變化，即沉管升沉方向與駁船具有相同的位移變化。在t+Δt 時刻，僅考慮駁船升沉運動影響時沉管瞬時位置ξzj（t+Δt）為

駁船的升沉運動使得沉管位置由ξj（t）變成ξzj（t+Δt），沉管受到的波浪力也發(fā)生變化。數值求解對應于沉管位置ξzj（t+Δt）的速度勢φS，帶入波浪力公式（12），可得到沉管在位置ξzj（t+Δt）所受到的波浪力Fk′。然后分別用ξzj（t+Δt）和Fk′代替公式（16）中的ξj（t）和Fk（t），即可得到t+Δt 時刻考慮波浪作用與駁船升沉運動共同影響的沉管位置ξj（t+Δt）。

3 算例計算及分析

數值計算中取沉管尺寸的長×寬×高為100 m×15 m×10 m，沉放深度h 為8 m，水深H 為40 m。入射波為正向規(guī)則波，周期T 為5 s，波高為1.0 m。纜繩為鋼纜，直徑為0.3 m，纜繩彈性系數為2.75×106kg/cm2。駁船升沉運動的振幅與相位取值見表1。

表1 計算工況Tab.1 Calculated condition

首先對忽略駁船運動時的沉管運動響應進行了計算。圖2 給出了忽略駁船升沉運動情形下沉管升沉運動時間過程線與入射波浪的時間過程線的比較，時間段取100～110 s（20～22 T）。圖中實線為忽略駁船運動情形下沉管升沉運動時間過程線，虛線為波浪的時間過程線。由圖2 可見，忽略駁船運動時的沉管升沉運動幅值約為0.1 m，沉管升沉運動過程線與入射波浪過程線相比，兩者之間存在較大的相位差，相位差約為π/2。

圖2 忽略駁船運動情形的沉管升沉運動時間過程線（H=1.0 m，T=5 s）Fig.2 Time history of waves and heave motion of tunnel element without barge motion（H=1.0 m，T=5 s）

圖3 不同相位情形下駁船升沉運動過程線（b=0.3 m）Fig.3 Time history of barge heave motion under different phases（b=0.3 m）

參照表1 給出的計算工況，對同時考慮波浪作用與駁船升沉運動影響的沉管運動響應進行了數值計算與分析。圖3 是駁船升沉運動振幅為0.3 m、駁船運動與波浪相位差分別取-π/2，-π/5，0，+π/5 和+π/2 情形的駁船升沉運動歷時過程線，時間段為100～110 s（20～22 T））。由圖3 與圖2 可見，當駁船升沉運動與入射波浪之間的相位θ=-π/2 時，駁船升沉運動與忽略駁船運動時的沉管升沉運動之間的相位差為0；當駁船升沉運動與入射波浪之間的相位θ=π/2 時，駁船升沉運動與忽略駁船運動時的沉管升沉運動之間的相位差為π。

圖4 是駁船升沉運動振幅為0.3 m、相位分別為-π/2，0 和+π/2 情形，考慮波浪作用和駁船升沉運動共同影響的沉管運動響應分量與忽略駁船運動時的沉管運動響應分量的比較。圖中heave 表示沉管升沉方向運動，sway 表示沉管橫蕩方向運動，roll 表示沉管縱搖方向運動。

由圖4 可見，當駁船做升沉運動時的沉管運動過程線與忽略駁船運動時的沉管運動過程線（圖中實線）相比，考慮駁船升沉運動時的沉管運動響應與忽略駁船運動時的沉管運動響應相比，在運動幅值與相位上都發(fā)生較大的變化。在sway 方向和roll 方向，考慮駁船升沉運動時的沉管運動響應的歷時過程線沒有變化，表明駁船升沉運動對沉管在sway 和roll 方向的運動影響很小。

圖4 駁船做升沉運動與忽略駁船運動情形的沉管運動響應（b=0.3 m）Fig.4 Motion response of tunnel element with and without barge heave motion（b=0.3 m）

為進一步考察駁船升沉運動相位差的影響，圖5 給出了駁船升沉運動幅值為0.3 m、駁船運動與波浪之間的相位差分別取-π/2，-π/5，+π/5 和+π/2 情形的沉管升沉運動過程線。圖中實線為考慮波浪作用與駁船升沉運動共同影響的沉管升沉運動時間過程線，虛線為忽略駁船升沉運動僅受波浪影響的沉管升沉運動時間過程線。為便于討論，駁船升沉運動的時間過程線也在圖5 中給出（圖中的點劃線）。由圖5 可見，當駁船升沉運動的相位為(-π/2)時，沉管升沉方向的運動幅值最大約為0.4 m，此情形下駁船升沉方向運動的相位與忽略駁船運動時的沉管升沉方向運動的相位一致，沉管升沉方向的運動幅值為駁船升沉運動幅值（0.3 m）與忽略駁船運動時的沉管升沉方向運動幅值（0.1 m）之和；當駁船升沉運動的相位為π/2 時，沉管運動幅值最小約為0.2 m，此情形下駁船升沉方向運動的相位與忽略駁船運動時的沉管升沉方向運動相位剛好相反，沉管升沉方向的運動幅值為駁船升沉幅值與忽略駁船運動時的沉管升沉方向運動幅值之差。

圖5 駁船升沉運動相位不同情況下沉管升沉運動過程線（b=0.3 m）Fig.5 Time history of tunnel element in heave motion under barge heave motion with different phases（b=0.3 m）

圖6 是駁船升沉運動振幅b=0.1～0.5 m、相位θ=-π/2～π/2 情形對應的沉管升沉運動幅值（H=1.0 m，T=5 s）。

由圖6 可見，當駁船升沉運動的相位θ 一定時，沉管運動幅值隨著駁船升沉運動的振幅增加而相應增加。駁船升沉運動的相位θ=-π/2 時，對不同振幅的駁船升沉運動，沉管升沉方向的運動幅值最大值為駁船升沉運動幅值與忽略駁船運動時沉管升沉運動幅值之和，發(fā)生在駁船升沉運動相位與忽略駁船運動時沉管升沉方向運動的相位一致的情形。駁船升沉運動的相位θ=π/2 時，對不同振幅的駁船升沉運動，沉管升沉方向的運動幅值最小值為駁船升沉運動幅值與忽略駁船運動時沉管升沉運動幅值之差，發(fā)生在駁船升沉運動相位與忽略駁船運動時沉管升沉方向運動的相位剛好相反的情形。

圖6 駁船做不同升沉運動情況下的沉管升沉運動幅值（H=1.0 m，T=5 s）Fig.6 Movement amplitude of tunnel element in heave motion under barge with heave motion（H=1.0 m，T=5 s）

4 結論

本文考慮駁船升沉運動與波浪的共同作用，建立了時域求解沉管運動的數學模型，分析駁船做不同振幅與不同相位的升沉運動對沉管運動特性的影響規(guī)律，得到如下主要結論：

（1）與忽略駁船升沉運動時的沉管運動相比，駁船升沉運動振幅與相位對沉管升沉運動有較大的影響。當駁船升沉運動的相位θ 一定時，沉管運動幅值隨著駁船升沉運動的振幅增加而相應增加。

（2）沉管升沉方向的運動幅值最大值為駁船升沉運動幅值與忽略駁船運動時沉管升沉運動幅值之和，發(fā)生在駁船升沉運動相位與忽略駁船運動時沉管升沉方向運動的相位一致的情形。沉管升沉方向的運動幅值最小值為駁船升沉運動幅值與忽略駁船運動時沉管升沉運動幅值之差，發(fā)生在駁船升沉運動相位與忽略駁船運動時沉管升沉方向運動的相位剛好相反的情形。

（3）在實際工程中采用駁船沉放法進行沉管管段沉放施工作業(yè)中，可通過改變駁船升沉運動的相位來減小沉管管段升沉運動。本文的研究工作是在假定駁船做已知的升沉運動的情形下進行的，實際上駁船運動和沉管運動需要在波浪場內耦合求解，該問題有待于進一步的研究。

［1］ZHAN D X，WANG X Q. Experiments of hydrodynamics and stability of immersed tube tunnel on transportation and immersing［J］. Journal of Hydrodynamics：Ser. B，2001（2）:121-126.

［2］Jensen O P，Olsen T H，Kim C W，et al. Construction of immersed tunnel in off-shore wave conditions Busan-Geoje project South Korea［J］. Tunneling and Underground Space Technology，2006，21:333.

［3］陳智杰.波浪作用下沉管管段沉放運動的試驗與數值研究［D］.大連:大連理工大學，2009.

［4］Hess J L，Smith A M O. Calculation of non-lifting potential flow about arbitrary three-dimensional bodies［J］. Journal of Ship Research，1964，8（2）:22-44.

［5］GARRISON C J. Hydrodynamics loading of large offshore structures: Three-dimensional source distribution method［C］//ZIENKIEWICZ O C，LEWIS R W，STAGG K G W. Numerical Methods in Offshore Engineering. Chichester：wiley，1978:87-140.

［6］Liu Y H，Kim C H，Lu X S. Comparison of higher-order boundary element and constant panel methods for hydrodynamic loading［J］. International Journal Offshore and Polar Engineering，1991，1（1）:8-17.

［7］Teng B，Eatock Taylor R. New higher-order boundary element methods for wave diffraction/radiation［J］. Applied Ocean Research，1995，17（2）:71-77.

［8］Wehausen J V，Laitone E V. Surface waves［M］. Berlin: Springer-Verlag，1960：446-778.

［9］Stoker J J. Water Waves，Pure and applied mathematics［M］. New York: Inter science Publishers，Inc.，1957.