水下爆炸及隨機外力下破損艦船的傾覆概率分析

張 婧，施興華

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

1 引 言

破損艦船的生命力問題因其具有重要的軍事戰(zhàn)略意義而受到國內(nèi)外許多學者的關注。作戰(zhàn)環(huán)境中的水面艦船，除了受到風浪等隨機擾動外，還受到其他非隨機的不規(guī)則擾動。艦船受到武器攻擊發(fā)生破損后，進水產(chǎn)生的固定橫傾角將大大降低其在海上的戰(zhàn)斗力與生命力，在下次攻擊武器的爆炸載荷以及風浪載荷聯(lián)合作用下，穩(wěn)性降低，可能造成其整體傾覆。從歷次世界大戰(zhàn)來看，具有強裝甲艦船強度的破壞概率遠小于整船的傾覆概率，確定該隨機事件的發(fā)生概率具有重要的意義。

艦船傾覆因涉及到外載荷的隨機性和大幅度橫搖的強非線性，而使問題十分復雜。波浪載荷和風載荷均具有隨機性，這就導致艦船傾覆是具有一定概率的隨機事件。確定該隨機事件的發(fā)生概率，已成為急需解決的問題。許多學者在這方面做了大量的研究[1-2]。但目前尚未見到關于水下爆炸載荷和隨機風浪聯(lián)合作用下破損艦船傾覆概率研究的文獻資料。

本文推導了破損艦船遭受水下爆炸及有色噪聲時的橫搖微分方程?？紤]到本文主要研究爆炸載荷對破損艦船橫搖的影響以及高維FPK方程求解繁瑣，故將有色噪聲過程簡化為白噪聲過程，以Gauss-Legendre路徑積分為基礎，給出橫搖角的概率密度函數(shù)隨時間演變的計算方法，按照現(xiàn)有的傾覆準則給出預報一定裝藥量下艦船傾覆概率的表達式，并進行了實例分析，對作戰(zhàn)環(huán)境中破損艦船傾覆研究進行有益的探索和初步嘗試。

2 作戰(zhàn)條件下破損艦船的環(huán)境載荷

2.1 水下爆炸載荷產(chǎn)生的傾側力矩

沖擊波隨時間衰減規(guī)律的經(jīng)驗公式為[3]

式中：θ—指數(shù)衰減的時間常數(shù)，對于TNT炸藥而言，θ=（W1/3/R）-0.24W1/3×10-4；pm—壓力峰值[4]，pm=53.3（W1/3/R）1.13，W—裝藥量，R—測點到爆炸中心的距離。

將作用于艦船水下舷側的爆炸沖擊波載荷視為一個突加傾側力矩

式中：A1—沖擊波作用面積，取塑性區(qū)域的面積，A1=π [（10～2 0） R0]2，R0—炸藥初始半徑；Me（t）—集中力 P（t）豎向偏心引起的傾側力矩；Δz1—集中力P（t）作用點離艦船重心的豎向距離。

將（1）式代入（2）式，得到

2.2 隨機風浪擾動力矩

風擾動力矩可以看作是平均風傾力矩Ma與脈動風傾力矩Md（t）兩部分組成，即

式中：ρ—空氣密度；A—艦船結構水線以上部分側投影面積；Δz—艦船結構水線以上部分側投影面積的形心至水壓力作用點的距離；U（t）—艦船的航速；υd（z,t）—短周期內(nèi)脈動風速，其譜密度Sυd（ω）由文獻[5]確定；υa（z）—平均風速[6]；Cm—風壓傾側力矩系數(shù)[7]，

式中：L—艦船結構總長；B—艦船橫剖面寬度；z—艦船結構水線以上部分側投影面積的形心至水線的距離。

（6）式充分反映了船型主要因素對風壓的影響，其變化范圍在0.955～1.418之間，如果Cm取為定值，勢必對有些船舶穩(wěn)性要求過高，而對有些船舶又會顯得不足。因此，根據(jù)船型選取不同的Cm值是較符合客觀實際的。

波浪擾動力矩為

航速U（t）與給定z處的平均風速都是常數(shù)，故由（4）式可知，Ma也是一常數(shù)。脈動風的強度隨時間而隨機變化，是典型的隨機過程。 由（5）式可知，Md（t）是 υd（z,t）的線性函數(shù)，所以可認為Md（t）是平穩(wěn)正態(tài)隨機過程。同理，Mw（t）也是平穩(wěn)正態(tài)隨機過程。

3 破損艦船橫搖運動微分方程

由于海浪及風力均為隨機過程，即使爆炸載荷是確定性載荷，但其作用時間短暫，所以作戰(zhàn)環(huán)境中艦船的運動也具有隨機性。由于爆炸載荷持續(xù)時間一般只有幾秒左右，平均風傾力矩作用時間一般在10min以上，所以在爆炸載荷作用過程中，平均風傾力矩始終作用在艦船上，可以把爆炸載荷產(chǎn)生的短時力矩及平均風傾力矩看作定常傾側力矩。則破損艦船大幅橫搖運動的隨機微分方程為

初始條件為

式中：φ01—遭多次攻擊后艦船破損進水產(chǎn)生的初始橫傾角。這里只考慮艦船的線性橫搖阻尼系數(shù)，則

式中：BL—線性阻尼系數(shù)。

采用線性項加三次項作為回復力矩模型,能較好地擬合靜穩(wěn)性曲線[10]，則

式中：Δ—排水量；GZ（ φ（t ））—回復力臂；C3—三次回復力矩系數(shù)，取為負值。

將（8）式兩邊同時除以 Jφφ+ΔJφφ，可得

式中：ωφ—為橫搖固有頻率，

隨機橫搖力矩Y（t）與Z（t）的統(tǒng)計性質(zhì)可分別由Mw（t）與Md（t）的統(tǒng)計性質(zhì)確定，自相關函數(shù)分別為 RY（τ）和RZ（τ），功率譜密度分別為

4 非線性隨機橫搖響應的時域分析

4.1 濾波系統(tǒng)的確定

在時域內(nèi)考察一個隨機過程，Markov過程理論是一個重要方法。構造兩個線性濾波器，其輸入為白噪聲過程，輸出則為滿足給定統(tǒng)計特征值的隨機過程。為模擬實際風浪，取形狀濾波器分別為

式中：α1，β1，γ1，β2，γ2—控制濾波函數(shù)特性的參數(shù)；Ni（t） （i=1,2 ）—單位白噪聲過程，N1（t）與N2（t）相互獨立，Ni（t）=dWi（t）/dt，Wi（t）—單位Wiener過程。

而由（14）式和（15）式確定的譜密度分別為：

則令（12）式與（16）式的峰值、卓越頻率、譜面積相等，就可得到濾波參數(shù) α1，β1，γ1，同理可得 β2，γ2。

聯(lián)立（11）、（14）式與（15）式，得到一個新的擴維運動微分方程組為

4.2 伊藤隨機微分方程

將微分方程組（18）改寫為狀態(tài)方程組的形式

這樣就建立了求解破損艦船在作戰(zhàn)環(huán)境下橫搖的隨機微分方程組。

4.3 相應的簡化方程

實際隨機風浪必須處理為有色噪聲，理論上是可以通過增加濾波系統(tǒng)并利用Markov過程理論來求解這一復雜的問題，但這將出現(xiàn)高維的FPK方程，至今超過三維的FPK方程的計算量已經(jīng)非常耗費計算機資源?？紤]到本文主要是研究爆炸載荷對艦船傾覆概率的影響，由于爆炸載荷極短時間內(nèi)就達到最大值，此時平均風傾力矩與其相比非常小，所以只考慮脈動風的作用，且爆炸載荷強度遠遠大于隨機風浪擾動的強度，在不影響問題結論的基礎上，可將隨機風浪的有色噪聲簡化為白噪聲過程。即將（19）式簡化為

式中：X1（t）=φ（t ）；X2（t）=（t）；γ—隨機擾動力矩系數(shù)；N（t）—單位白噪聲過程。

初始條件為

通常上述方程組的求解是通過其相應的偏微分方程進行的， 即轉(zhuǎn)移概率密度）滿足FPK方程

4.4 路徑積分法求解方程

路徑積分的基本思想就是在空間和時間上分別離散化，以路徑和代替積分，即通過連接短時轉(zhuǎn)移概率密度形成全局轉(zhuǎn)移概率密度，得到狀態(tài)向量的聯(lián)合概率密度函數(shù)。設（t）是n維狀態(tài)方程，其演化概率密度為

路徑積分是在縮減的狀態(tài)空間Rs內(nèi)積分，Rs之外的區(qū)域轉(zhuǎn)移概率密度充分小，可忽略不計，即

Yu[11]等人將（23）式按照Gauss-Legendre積分來離散化，得到基于Gauss-Legendre公式的路徑積分法。已知第（i-1）時刻的每個高斯積分點上的概率密度及相應轉(zhuǎn)移概率密度時，借助離散化的概率密度表達式，可獲得第i時刻任意點的概率密度，因此僅需計算第i時刻相應高斯積分點的轉(zhuǎn)移概率密度。這樣可大大減少計算量，最終可得到適合編程實現(xiàn)的表達式。對于一維的情況，即

式中：K—子區(qū)間數(shù)；Lk—第k子區(qū)間的高斯積分點數(shù)；δk—第k子區(qū)間的長度；xkl—高斯積分點；ckl—相應的權數(shù)。

一般假設短時轉(zhuǎn)移概率密度是近似高斯分布的。Sun和Hsu[12]提出利用矩方程導出短時轉(zhuǎn)移概率密度的一階矩和二階矩，由于非線性隨機系統(tǒng)的矩方程一般是無窮層階次、非封閉的，可利用高斯截斷法獲得封閉的矩方程。因此，短時概率密度可寫為

式中：E1（ti），E2（ti）—ti時刻X的一階、二階原點距，可通過封閉的矩方程組[13]解出；σ2（ti）=E2（ti）-[E1（ti）]2。

5 破損艦船的傾覆概率分析

艦船結構傾覆是由于橫搖角過大引起的，當橫搖角超過某一閾值，就將一直增大而不會回到平衡位置。艦船的穩(wěn)性消失角是這一閾值的最好近似，因此將艦船結構傾覆定義為其橫搖角大于穩(wěn)性消失角或傾覆角，即

處于正浮狀態(tài)下的艦船，當回復力矩為零時，其橫搖角定義為穩(wěn)性消失角或傾覆角φv，即

得到方程的正根與負根為：

關于艦船傾覆問題的討論中，傾覆概率的定義不是很明確，艦船的傾覆是一個不可再現(xiàn)的事件，因此定義在某一固定時刻t，艦船結構的傾覆概率為

式中： f（ φ,t ）—橫搖角邊緣概率密度分布。

由前面的分析可知， f（ φ,t ）與時間及裝藥量密切相關，所以根據(jù)（26）式可預報任意裝藥量下任意時刻破損艦船的傾覆概率。

6 計算實例及分析

采用fortran90語言自編程序，對某艦船遭受爆炸載荷、白噪聲風浪及平均風傾力矩聯(lián)合作用下橫搖傾覆進行了分析計算。爆炸沖擊載荷的參數(shù)見表1。與隨機擾動及艦船有關的參數(shù)如表2所示。

表1 與水下爆炸載荷有關的參數(shù)Tab.1 The parameters of underwater explosions

（1）破損艦船再次受到攻擊時橫搖角概率分布隨時間的變化

破損艦船在裝藥量W分別為0kg，600kg，900kg時橫搖角的邊緣概率密度隨時間的演變?nèi)鐖D1-3所示。

表2 與隨機擾動及艦船有關的參數(shù)Tab.2 The parameters of random disturbance and warship

圖1 W=0kg時，t分別為0s，0.1s，1.5s時橫搖角的邊緣概率密度Fig.1 The marginal probability density of rolling angle at 0s,0.1s and 1.5s(W=0kg)

圖2 W=600kg時，t分別為0s，0.1s，1.5s時橫搖角的邊緣概率密度Fig.2 The marginal probability density of rolling angle at 0s,0.1s and 1.5s(W=600kg)

圖3 W=900kg時，t分別為0s，0.1s，1.5s時橫搖角的邊緣概率密度Fig.3 The marginal probability density of rolling angle at 0s,0.1s and 1.5s(W=900kg)

從圖1-3可看出，具有定常側傾力矩的艦船在隨機外力作用下橫搖角的概率分布發(fā)生漂移和擴散，并且不再圍繞正浮位置分布，而是圍繞定常側傾力矩造成的橫傾角分布。同一裝藥量下，由于隨機外力的影響，峰值隨時間的增長而降低，但橫搖角在較大的范圍內(nèi)變化。除了t=0s時刻橫搖角的邊緣概率密度分布不受裝藥量的影響外，其它時刻的峰值隨裝藥量的增加而降低，由定常傾側力矩造成的傾角增加。

在裝藥量W=0kg時，各分布圖所圍繞的固定橫傾角不變。其它裝藥量的情況時，同一裝藥量下，雖然隨著時間的增長，爆炸載荷產(chǎn)生的力矩減小，導致固定橫傾角呈減小的趨勢，但是前一時刻產(chǎn)生的橫傾角對后一時刻橫搖角的概率分布有累積影響，而且由于隨機外力的影響，橫搖角的分布逐漸擴散。

（2）裝藥量對破損艦船傾覆概率的影響

由以上分析可知，裝藥量對橫搖角的概率分布有明顯的影響，即分布的峰值降低，橫傾角增大。而且裝藥量的增加即固定橫傾力矩的增加將明顯增加傾覆的可能性。按（26）式計算不同裝藥量不同時刻破損艦船的傾覆概率，結果如表3所列。

表3 不同裝藥量不同時刻破損艦船的傾覆概率pfTab.3 The capsizing probability pfof damaged warship under different blasting charge at different time

從表3可以看出，對已破損艦船再次受到武器攻擊時，傾覆概率明顯增加。但由于爆炸載荷作用時間很短，故其后傾覆概率有所減小。

7 結 論

本文對水下爆炸載荷和隨機風浪聯(lián)合作用下破損艦船的橫搖運動及其傾覆概率進行了分析計算。由上述理論和實例分析，可得出如下結論：

（1）考慮到艦船的結構形式比較特殊，改變傳統(tǒng)的處理方法，將風壓傾側力矩系數(shù)Cm視為變化量，計入船型因素對風擾動力矩的影響，同時考慮了航速對風速的影響，即風以相對風速的形式作用于艦船，給出了計入航速影響的風擾動力矩計算公式，更符合實際情況。

（2）算例表明，水下爆炸載荷對橫搖角的概率分布有明顯的影響，即分布的峰值降低，橫傾角增大。而且裝藥量的增加即固定橫傾力矩的增加將明顯增加傾覆的可能性。已破損艦船再次受到武器攻擊時，傾覆概率明顯增加。

（3）本文給出了在水下爆炸及隨機風浪作用下破損艦船大幅橫搖運動微分方程，研究了隨機外力作用下船舶的一維運動。雖然這是與傾覆關系最密切的運動模態(tài)，但畢竟尚未給出船舶運動的全面描述，和其它運動模態(tài)耦合后的艦船橫搖的研究將作為今后進一步的研究方向。

（4）對于同時遭受爆炸載荷、隨機風力及海浪的破損艦船，本文引入兩個濾波器，將隨機動力學領域理論上較為成熟的路徑積分法推廣應用于求解其傾覆概率，即將隨機微分方程升到五維，通過求解相應的FPK方程得到橫搖角的概率密度函數(shù)，由于FPK方程的系數(shù)與風浪載荷及爆炸載荷密切相關，這樣便可以估計風浪載荷及爆炸載荷對艦船傾覆概率的影響。但目前，五維FPK方程還無法求解，這里給出了分析計算傾覆概率的思路，高維FPK方程的數(shù)值求解有待于以后進一步的研究。

[1]Barone G,Navarra G,Pirrotta A.Probabilistic response of linear structures equipped with nonlinear damper devices(PIS method)[J].Probab Eng Mech,2008,23(2-3):125-133.

[2]Cottone G,Di Paola M,Ibrahim R.Stochasic ship roll motion via path interal method[J].Inter J Nav Oc Engng.,2010,2:119-126.

[3]張 挺.爆炸沖擊波測量技術（電測法）[M].北京:國防工業(yè)出版社,1984.

[4]Rajendran R,Lee J M.Blast loaded plates[J].Marine Structures,2009,22:99-127.

[5]Davenport A G.The spectrum of horizontal gustiness near the ground in high winds[J].J Royal Meteorol Soc.,1961,87:194-211.

[6]張相庭.結構風壓和風振計算[M].上海:同濟大學出版社,1985.

[7]湯忠谷,韓久瑞.海船風壓試驗研究[J].中國造船,1981(2):31-38.

[8]李積德.船舶耐波性[M].哈爾濱:哈爾濱船舶工程學院出版社,1992.

[9]俞聿修.隨機波浪及其工程應用[M].大連:大連理工大學出版社,2003.

[10]Nayfeh A H,Khdeir A A.Nonlinear rolling of ship’s in regular beam seas[J].ISP,1986,33:40-49.

[11]Yu J S,Cai G Q,Lin Y K.A new path integration procedure based on Gauss-Legendre scheme[J].International Journal of Non-Linear Mechanics,1997,32(4):759-768.

[12]Sun J Q,Hsu C S.The generalized cell mapping method in nonlinear random vibration based upon short-time Gaussian approximation[J].Journal of Applied Mechanics,1990,57:1018-1025.

[13]方 同.工程隨機振動[M].北京:國防工業(yè)出版社,1995.