李林斌,李洛東,劉日明
(中國船級社,北京 100007)
準確預報船舶及海工結構在水動力作用下的長期極值響應是結構設計的基礎,如船舶25年重現期的波浪彎矩或半潛式平臺100 年重現期的氣隙等。長期響應預報需要考慮兩個因素,一是結構在穩定海況條件下響應的短期分布,二是海況條件的長期分布,進而將每個海況條件對結構長期極值響應的貢獻率統計在內。對于非線性響應的情況(如波浪砰擊顫振誘導的船舶垂向彎矩),需要針對每個海況條件,采用非線性時域模擬和統計分析方法,來獲得結構的響應分布,因此非常耗費計算資源,不便于結構的設計優化和敏感性分析。
采用合理的替代方法,避免在整個海況域內進行非線性響應分析,是船舶和海洋工程領域持續研究的課題。其中,環境包絡線法在海洋工程領域有較為普遍的應用[1],這是一種基于短期極值來推算長期極值的方法,即在特定重現期環境包絡線設計點海況的短期最大響應基礎上,乘以一系數(一般為1.1~1.3),或取短期最大響應分布的更高分位值(一般為75%至90%分位),作為相同重現期的長期極值響應,但該系數或分位值取決于結構響應的短期變異程度以及海況條件的長期變異程度,需要預先研究確定;另外一種方法是根據各海況對線性響應長期極值的貢獻度,篩選出一定數量的海況,再針對這些特定海況進行更為復雜的非線性數值模擬和統計分析[2],該方法成立的前提是海況的非線性和線性響應貢獻度排序大致相當,但這種假設并非普遍成立,同時并非所有的結構線性響應都存在且容易獲得。近年來,采用機器學習技術的回歸算法得到越來越多的重視和應用,該方法一般采用平均取樣[3]或拉丁超立方體抽樣技術[4],選取一定數量的海況進行水動力分析,獲得結構響應的分布參數,用以訓練代理模型,并基于代理模型來預測其他海況對應的響應分布參數。由于需要保證代理模型在整個輸入域內預測結果的準確度,訓練集樣本數往往很大,導致非線性響應分析的計算量依然巨大。總之,這些替代性方法的應用都存在一定的局限性。
本文提出一種新方法進行非線性水動力響應的長期預報。該方法采用機器主動學習技術,通過機器學習模塊和水動力計算模塊的交互,以循環迭代的方式,貫序篩選對極值響應貢獻最大的海況進行非線性時域分析,并將分析結果追加到訓練集中,持續訓練機器學習模型,進而實現預測值向真實值的快速收斂。
在后續的章節中,作者將簡要介紹采用傳統方法進行非線性水動力響應長期預報的過程;闡述如何利用高斯過程回歸(GPR,Gaussian Process Regression)和貫序取樣(Sequential Sampling)技術來建立機器主動學習模型;并以兩型大型集裝箱船為算例,采用機器學習模型對彈性體砰擊顫振彎矩和剛體非線性波浪彎矩進行預測;通過與真實值的比較,驗證模型的高效性和準確性,展示其在工程上的應用價值;最后給出幾點研究結論和建議。
船舶和海工結構在水動力作用下的長期極值響應,無論是線性的還是非線性的,都需要考慮結構在穩定海況條件下的響應短期分布,以及海況條件的長期分布,即
式中:Q(y)為響應Y的長期超越概率,目標超越概率取決于重現期的要求;f(x)為海況x的長期分布,如波高-周期(Hs-Tz)聯合概率分布函數或波浪散布圖;g(y|x)為響應Y在海況x下的短期分布(超越概率)。其中,響應的短期分布g(y|x)是整個過程的核心。對于線性響應,問題相對簡單,可以通過傳遞函數和譜分析的方法得到;對于非線性響應,就需要針對每個海況,通過非線性時域水動力計算和統計分析獲得,需要在整個海況域內求解,因此是一個非常耗時的過程。下面以船舶波浪砰擊顫振為例進一步說明。
大型船舶,特別是大型集裝箱船,由于首尾部船體大幅外飄,當以較高航速在中高海況航行時,大幅的船波相對運動會使船體受到波浪的強烈沖擊,船體結構由于砰擊載荷作用而發生瞬時高頻振動,產生所謂的“顫振”現象,具有強非線性特征,對總體波浪彎矩影響顯著,需要在船體梁極限強度分析時充分考慮。根據中國船級社《船體結構波激振動和砰擊顫振直接計算評估指南》,應采用時域水彈計算和統計分析方法進行預報,包括彈性體砰擊顫振誘導彎矩和剛體非線性波浪誘導彎矩。
對無限航區船舶,一般采用IACS REC No.34 推薦的北大西洋海浪散布圖[5],其海況數為197 個(即:197個非零Hs-Tz參數對)。浪向角Hd的選取原則是在0°至360°之間,取間隔不大于30°,各浪向角的發生概率相等。考慮到對稱性的因素,浪向角合計為7 個,如圖1 所示。這樣在考慮浪向后的海況總數為1379個,對如此多的海況進行逐一的非線性時域分析,對計算資源的耗費是可想而知的。

圖1 浪向角Hd及編號(括號內的數字為Hd編號)Fig.1 Wave heading Hd and its numbering(the digit in brackets is the Hd numbering)
在短期海況xji(第j海況第i浪向角)條件下進行非線性時域計算,對得到的垂向波浪彎矩峰值y進行Weibull分布擬合,其超越概率為
式中,βji為形狀參數,ηji為尺度參數。
垂向波浪彎矩長期極值超越概率為
式中,pj為第j海況的長期概率,pi為第i浪向的長期概率。Q( )y= 10-8所對應的響應值yex為目標極值響應,即25年重現期的垂向波浪彎矩,本文中即為彈性體的砰擊顫振彎矩或剛體的非線性波浪彎矩。
海況xji對長期極值響應的貢獻度Cji是一個很重要的概念,在后續的分析中會反復用到,其定義為
圖2 顯示了分析的總體過程。對Cji的進一步分析,不難發現長期極值響應來自于極少數海況的貢獻。以13 500 TEU顫振中垂彎矩為例(具體見3.1節),貢獻度前30的海況,其累計貢獻度為98.3%,而這些海況僅占海況總數(1379)的0.008%,均為有義波高Hs超過12.5 m 的高海況,浪向都是迎浪或艏斜浪(1號和2號浪向角),總體情況如圖3所示。

圖2 時域非線性波浪載荷計算和統計分析流程Fig.2 Time domain nonlinear wave load calculation and statistical analysis process

圖3 海況對長期極值響應的貢獻度(圖中圓的直徑與貢獻度成正比,顏色對應浪向角)Fig.3 Contributions of sea states to the long-term extreme response(the diameter of the circles is proportional to the contributions,and the colors correspond to the wave headings)
如果能夠建立一種模型,以智能和貫序方式識別出這些海況,并傳遞給水動力分析模塊進行非線性時域分析,而不是分析全部海況,將大大減少計算資源的耗費,提高結構設計和分析效率,這正是本文研究的目的和出發點。
高斯過程回歸(GPR)是近年發展起來的一種機器學習回歸方法[6],有著嚴格的統計學習理論基礎,對處理高維數、小樣本、非線性等復雜的問題具有很好的適應性,且泛化能力強,有著廣泛的工程應用[3,7]。GPR從函數空間角度出發,定義一個高斯過程來描述函數分布,直接在函數空間進行貝葉斯推理。GPR 是任意有限個隨機變量均具有聯合高斯分布的集合,其性質完全由均值函數和協方差函數確定,即f(x)~GP(μ(x),k(x,x') ),其中
式中,x,x'∈Rd為任意隨機變量。通常會對數據作預處理,使其均值函數等于0。假設觀測過程受到加性噪聲的污染,則觀測值y可用如下模型描述:
式中,x為輸入向量,f為函數值,ε為加性噪聲污染。假設ε符合均值為0、方差為c2的正態分布,即ε~N(0,c2),則觀測值y的先驗分布為
觀測值y和預測值f*的聯合先驗分布為
式中:K(X,X)=Kn=(kij)為n×n階協方差矩陣,kij=k(xi,xj),為xi與xj的協方差;K(X,x*)=K(x*,X)T為測試點x*與訓練集的輸入X之間的n× 1 階協方差矩陣;k()x*,x*為測試點x*的方差;In為n維單位矩陣。
由此可以得到預測值f*的后驗分布為
為了訓練GPR模型,需要設定一個協方差函數,文中我們采用了高斯徑向基核函數[3]:
式中,δ2為核函數的方差參數,h為尺度參數。δ和h這兩個超參數以及c需要通過機器學習(即訓練)獲得。
GPR 是監督學習的一種,是根據經驗數據(訓練集)來學習輸入-輸出之間的映射關系,使得給定新的輸入便可得到相應的輸出值(預測值)。通常情況下,訓練集需要人為預先指定。相比人工神經網絡、隨機森林、支持向量機等其他機器學習方法,高斯過程回歸的獨到之處在于其預測值是一個高斯隨機變量,不但可以給出預測值的最佳估計(均值),而且可以給出變異程度(方差)。
正是GPR以上獨有的特點,為實現機器主動學習提供了可能,并應用于非線性動力系統的長期極值響應預報[8],即貫序找出對長期響應影響最大的輸入值,將其加入到訓練集,迭代更新GPR 模型,使其預測值在對長期響應貢獻度最大的輸入區域內變得更為準確,進而實現預測值向真實值的快速收斂。結合本文的研究目的,我們設計了一種算法,通過Python編程,開發了機器主動學習模塊,并通過與水動力分析模塊的交互,實現了對船舶結構非線性水動力響應的長期預報,具體步驟如下:
(1)指定初始訓練集的輸入集,即海況x=()Hs,Tz,Hd,進行非線性時域水動力計算和統計分析,得到訓練集的輸出集,即Weibull分布的形狀參數β和位置參數η;
(2)對GPR模型進行訓練,分別獲得超參數δβ、hβ和cβ以及δη、hη和cη;
(3)利用GPR模型預測其他海況xji對應的觀測值均值和,以及方差和;
Murphy指出,個體在實際展開工作活動的過程中,為了達到工作績效目標所表現出的一系列行為,即為績效;我國在積極展開工作績效研究的過程中,指出了學習過程、創新行為、公民氣候、技術核心即學習績效、創新績效、關系績效和任務績效四個工作績效結構[2]。
(5)計算每一個海況的“收獲函數”Aji=,找出Aji最大時對應的海況xn;
(6)將xn傳遞給水動力分析模塊,進行非線性時域計算和統計分析,獲得Weibull 分布的形狀參數βn和位置參數ηn;
(7)將xn,βn和ηn追加到訓練集中,重復(2)~(6),直到長期響應值收斂。
以上過程中,“下一個最優”海況的選擇應考慮和平衡兩個因素,一是需要距離訓練集足夠遠(也即預測值方差足夠大),二是對長期極值響應y^ex的貢獻度足夠大,而“收獲函數”實現了這一平衡。
本文以13 500 TEU 和20 000 TEU 兩型大型集裝箱船為算例進行研究,建立水動力模型并進行時域水彈計算和統計分析,包括彈性體砰擊顫振彎矩和剛體非線性波浪彎矩。這兩型集裝箱船的主要參數如表1所列。

表1 船舶主要參數Tab.1 Main parameters of the ships
首先采用傳統方法,按照1.2 節所述流程,完成了全部1379 個海況的水動力計算和統計分析,并獲得了長期極值響應yex,并將其作為“真實值”。
然后采用第2 章所述的機器學習模型進行長期極值響應預報。選擇波浪散布圖中最小和最大Hs,及其對應Tz出現概率最大的海況,對這些海況進行水動力計算和統計分析,獲得對應的分布參數,作為初始訓練集。
式(13)為13 500 TEU船型顫振中垂彎矩GPR模型的初始訓練集:
我們采用機器學習模型迭代100次,每次迭代都獲得一次長期極值響應的預報值y^ex;同時模型將根據“收獲函數”,自動和貫序識別出“下一個”最重要的海況,傳遞給水動力模塊進行時域分析。給出的誤差e如式(14)所示,以反映機器學習模型的準確性:
誤差分析的結果如圖4~7 所示。不難看出,機器學習模型的預測結果誤差隨著迭代次數的增加而快速降低和收斂。雖然由于結構響應的非線性、統計過程的不確定性,以及機器學習過程的隨機性等因素導致不同船型的收斂效率會有所差別,但從本文的算例看,模型迭代大約30~40 次后,其誤差均可以控制在1.5%以內。加上初始訓練集的14個海況,需要進行非線性水動力分析的海況總數不超過60個,相比傳統方法需要分析的全部1379個海況,計算效率提高了20多倍。

圖4 13 500 TEU彈性體砰擊顫振彎矩誤差Fig.4 Errors of the 13 500 TEU elastic hull whipping bending moments

圖5 13 500 TEU剛體非線性波浪彎矩誤差Fig.5 Errors of the 13 500 TEU rigid body nonlinear wave bending moments

圖6 20 000 TEU彈性體砰擊顫振彎矩誤差Fig.6 Errors of the 20 000 TEU elastic hull whipping bending moments

圖7 20 000 TEU剛體非線性波浪彎矩誤差Fig.7 Errors of the 20 000 TEU rigid body nonlinear wave bending moments
圖8 顯示了機器學習模型貫序取樣得到的前20 個海況的排序(13 500 TEU 彈性體顫振中拱彎矩)。與圖3 對比可以看出,機器學習模型很好地識別出了對長期極值響應貢獻度最大的海況區域,顯示出機器學習技術的強大功能。

圖8 貫序取樣前20個海況的分布,0表示初始訓練集(13 500 TEU顫振中拱彎矩)Fig.8 Scatter of the top 20 sea states by sequential sampling and 0 means the initial training dataset(the 13 500 TEU whipping hogging moment)
第2章所述的機器學習模型中,初始訓練集的輸入集(海況)需要人為預先設定,并對收斂效率產生影響。海況涉及到三個參數,即:有義波高Hs、跨零周期Tz和浪向角Hd。我們以13 500 TEU 顫振中垂彎矩為例,通過分別改變式(13)中的這三個海況參數,進行敏感性分析,結果匯總如圖9所示。

圖9 訓練集的敏感性分析(13 500 TEU彈性體顫振中垂彎矩)Fig.9 Sensitivity analysis of training data(the 13 500 TEU elastic hull whipping sagging moments)
(1)改變有義波高
一般情況下,船舶結構的水動力長期極值響應主要來自于大波高海況(如Hs超過10 m 的海況),圖3 也清晰反映了這一趨勢。因此,在確定初始訓練集的下限海況時,似乎應選擇更高的海況,以便提高收斂效率和計算精度,但事實并非如此。如我們選擇(Hs,Tz)=(9.5,10.5)作為下限海況,上限海況保持不變,可以看到模型的收斂效率明顯下降,原因在于GPR 是一種適于內插預測的機器學習方法,而外推預測時誤差往往會很大,即模型對Hs低于9.5 m海況的貢獻度預測產生很大誤差,進而導致收斂效率降低。
(2)改變波浪周期
選擇波浪散布圖中的最小波高和最小周期(Hs,Tz)=(0.5,3.5),以及最大波高和最大周期(Hs,Tz)=(16.5,15.5),作為初始訓練集海況,其收斂效率會有所下降,但不是很明顯。
(3)改變浪向角
只保留1/4/7三個浪向,將式(13)的海況訓練集樣本數縮減至6個。可以看出,誤差收斂效率有所降低,需要迭代更多次才能達到預期的誤差范圍。
根據以上分析,在確定初始訓練集海況時,我們應選擇波浪散布圖中最小和最大Hs及其對應Tz出現概率最大的海況,以及所有的浪向角,以保證迭代效率和預測精度。
船舶和海工結構在水動力作用下的響應往往具有非線性特征,對長期響應極值的準確預報通常需要在整個海況域上進行時域非線性數值模擬和統計分析,計算資源耗費巨大。本文采用高斯過程回歸和貫序取樣技術,設計了機器主動學習的算法,并通過編程開發出了相應的應用模型,實現了基于少數關鍵海況的分析結果來準確預測長期響應極值的目標。
本文以兩型大型集裝箱船為算例,采用機器學習模型,對彈性體砰擊顫振彎矩和剛體非線性波浪彎矩的長期響應極值進行了預報,結果表明其預測值隨著迭代次數增加而快速收斂,在可接受的誤差范圍內,較傳統方法提高計算效率20多倍,達到了工程化應用的程度。對貫序取樣結果的分析表明,模型可以很好地識別出對極值響應影響最大的輸入區域,顯示了機器學習技術的強大功能。
以上機器學習模型中沒有引入特定假設,因此具有廣泛的適用性,做稍許適應性改造就可以用于其他非線性響應的長期預報。如果將這一機器學習模塊嵌入到水動力分析軟件中,實現模型迭代結果和水動力計算結果的自動交互,將大幅提高軟件的智能化水平和計算效率。