基于Sobol’法的轉輪葉片動態響應靈敏度分析

李兆軍，張光政，劉福秀，丁 江

（廣西大學機械工程學院，廣西南寧 530004）

0 引言

水輪發電機組在運行過程中的工況復雜多變，致使水流常常處于非穩定狀態，從而使轉輪葉片在復雜的水流激勵作用下產生劇烈振動，影響水輪機的安全穩定運行［1］，因而實現轉輪葉片振動特性的有效調控使其振動特性滿足水輪機正常工作是十分必要的。目前，學者們對轉輪葉片振動特性控制方面的研究主要是集中在如何避免轉輪葉片產生水力共振等方面。例如，史廣泰等［2］采用ANSYS 平臺對水輪機轉輪進行了模態計算和振動分析，分析了機組產生共振的運行條件；李海亮等［3］采用有限元法研究了混流式水輪機轉輪振動特性以及各部件對轉輪固有頻率的影響規律；張新等［4］為了使轉輪的固有頻率避開外界的激振頻率，利用CFD 軟件對水輪機轉輪進行了有限元模態分析。但這些研究均未討論轉輪葉片振動特性對工況參數的敏感程度。

動態響應靈敏度能夠反映設計變量或者參數對結構振動指標的影響，許多學者開展了相關方面的研究，例如，孫維光等［5］基于Sobol’法推導了動力包機組振動對懸掛隔振參數的全局靈敏度；王長科等［6］對列車車體的扭轉振動進行了靈敏度分析，得到了對車體扭轉振動頻率影響較大的設計參數；劉柏希等［7］利用偏導數法分析了RV 減速器系統固有頻率對轉動慣量和剛度的影響。因而，積極開展水輪機轉輪葉片動態響應靈敏度方面的研究對于有效控制轉輪葉片振動特性具有十分重要的意義。然而，目前關于轉輪葉片動態響應靈敏度問題的研究鮮見文獻報道。

本文以水輪機轉輪葉片為研究對象，根據轉輪葉片的動力學模型，得到轉輪葉片的動態響應表達式，在此基礎上利用Sobol’法推導轉輪葉片動態響應對工況參數的靈敏度，揭示工況參數對轉輪葉片動態響應的影響，并通過實例進行分析。

1 轉輪葉片動力學模型

考慮到轉輪葉片的幾何形狀以及在工作過程中的受力情況，其往往會發生彈性變形，故將轉輪葉片視為彈性體，采用矩形平面殼單元進行轉輪葉片的模擬。利用有限元法，建立轉輪葉片的動力學方程為［8］：

式中：M為轉輪葉片整體質量矩陣，且M=Mt+Mp；Mt為單元的質量矩陣；Mp為附加質量矩陣，與水力參數有關，因而質量矩陣M不僅與轉輪葉片的結構參數有關，而且還和水力參數有關，為流固耦合項；K為轉輪葉片整體剛度矩陣；C為轉輪葉片阻尼矩陣；、、U分別為系統廣義位移向量、廣義速度向量、廣義加速度向量；Mtr為轉輪葉片的自激慣性力項，它不僅與轉輪葉片的運動參數有關，而且還與轉輪葉片的結構參數和水力參數有關，為流固耦合項；F1為卡門渦引起的壓力脈動，可由文獻［9］確定，F2為尾水管渦引起的壓力脈動，可由文獻［10］確定，F3為動靜干涉引起的壓力脈動，可由文獻［11］確定。

根據動力學方程，轉輪葉片頻率方程可表示為：

根據式（2），即可求出轉輪葉片的固有頻率。進一步，利用振型疊加法求解出轉輪葉片的振動響應：

式中：N為轉輪葉片自由度數；Φ為正則模態矩陣；η為相應的模態坐標列陣。

2 動態響應靈敏度分析

靈敏度分析可以反映輸入變量對輸出變量的影響程度，從數學上可表示為：若某一函數f（x）可導，則其一階靈敏度可表示為［12］：

式中：xi為所分析參數變量的集合；S又稱f（x）對xi的一階微分靈敏度。

這里采用Sobol’法來推導轉輪葉片動態響應對工況參數的靈敏度，Sobol’法是一種基于方差的蒙特卡洛法，其大致過程如下所示［13］：

首先，將轉輪葉片動態響應表示為U=f(X)，其中X=(x1，x2，…，xk)為影響轉輪葉片動態響應的k個工況參數所組成的參數向量，進而可建立輸入參數的空間域為Pk。

其次，將f（X）表示為2k個遞增項之和，即：

式中：f0為常數。其余子項對其所包含的變量的積分一定為零，即：

根據式（5）和（6），可以看出各項之間都是正交的，如果(i1，i2，…，is) ≠(j1，j2，…，jn)，有：

進一步有：

對式（5）取多重積分，即可得到各子項的表達式：

式中：X-i表示xi之外的其他工況參數，同理X-(ij)表示xi和xj之外的其他工況參數。

基于上述推導，可得到函數的總方差和偏方差的表達式為：

最后，各工況參數變量對轉輪葉片動態響應的靈敏度可表示為：

式中：Di1，i2，…，is和D分別為f（x）的偏方差和總方差；Si1為轉輪葉片動態響應對工況參數x1一階靈敏度，反映x1單獨變化對轉輪葉片動態響應的影響，同理可得其他轉輪葉片動態響應對其他工況參數變量的一階靈敏度。

3 實例分析

以國內某700 MW 大型混流式水輪機為例，主要計算參數如下：額定轉速n=107.14 r∕min；額定流量Q=554.52 m3∕s；額定負荷W=571 MW；轉輪最大外徑為D=8 051.8 mm；轉輪葉片數Z1=13，單個轉輪葉片的質量mt=7.7 × 103kg，轉輪葉片的材料密度為ρ=7.85×103kg∕m3，轉輪葉片的彈性模量為E=210 GPa，泊松比μ=0.3；導葉數Z2=26；轉輪葉片出水邊厚度δ=0.026 m；轉輪葉片出水口流角β=20°；斯特魯哈數sh=0.26；轉輪葉片尾部形狀影響系數k=2.4。考慮到水輪機工況參數較多，這里選擇一個易于控制的工況參數（轉速）進行分析，假設轉速的變化范圍為0.5n至1.2n，轉速每上升0.1n取一個轉速點進行分析，共8 個轉速點，記為ni（i=1，2，…，8），為了簡化計算，這里假設ni服從均勻分布，且取值范圍為［ni-0.02n，ni+0.02n］，如表1 所示。采用拉丁超立方采樣法對ni進行采樣，采樣次數設置為2 000次，并基于Sobol’法計算轉輪葉片動態響應對轉速的靈敏度，如圖1所示。

表1 轉速概率分布 r∕minTab.1 Probability distribution of speed

圖1 轉輪葉片動態響應對轉速的靈敏度曲線Fig.1 Sensitivity curve of blade dynamic response to speed

從圖1 中可以以看出，在轉速從n1變化至n8的過程中，轉輪葉片動態響應對轉速的靈敏度呈先增大后減小再增大的趨勢，當轉速變化至n6（0.6n）以及n8（1.2n）時，轉輪葉片動態響應對轉速的靈敏度出現了兩個峰值，說明在這兩個轉速附近轉輪葉片動態響應對轉速的變化是較為敏感的。

從上述分析可以看出，轉速對轉輪葉片動態響應的影響是隨轉速的大小而不斷變化的，如果要想實現在水輪機運行過程中對轉輪葉片振動的控制，應關注靈敏度較大的兩個位置，即轉速為n6和n8。

為了驗證仿真結果的正確性，進行了水輪機變轉速實驗，首先進行測點及裝置布置，即在轉輪葉片上布置應變片，并將應變信號采集儀安裝在主軸軸頂；然后讓水輪機分別在每個轉速下運行并采集動應力信號；最后進行數據處理并計算動應力通頻值，通過試驗采集得到轉輪葉片出水邊某一測點的動應力峰峰值隨轉速的變化曲線，如圖2所示。

圖2 轉輪葉片動應力隨轉速的變化曲線Fig.2 Variation curve of dynamic stress with rotating speed

從圖2 中可以看出，隨著轉速∕額定轉速的增大，轉輪葉片的動應力響應也隨之增大，當轉速比在［0.8，1］這一區間內，轉輪葉片動應力響應上升較為緩慢，處于一個較為穩定的狀態，總的來說，轉輪葉片動應力響應隨機組的轉速的增加而單調上升。

根據圖2 所示試驗所測得的轉輪葉片動應力隨轉速的變化曲線，計算其每個轉速比處的斜率，并進行歸一化處理，得到轉輪葉片動應力對轉速的靈敏度曲線，并將仿真得到曲線與其進行對比，如圖3所示。

從圖3 中可以看出，通過實測數據繪制成的靈敏度曲線和仿真所得到的靈敏度曲線基本一致，其峰值都出現在轉速為n6和n8處，只是靈敏度的大小有所差異。此外，該靈敏度曲線也與圖2 中轉輪葉片動應力隨轉速的變化趨勢較為符合，故可以說明本文通過仿真計算得到的結果是較為準確的。根據靈敏度分析的結果，可得到對轉輪葉片動態響應變化最敏感的轉速范圍，從而為轉輪葉片的振動控制提供依據。

圖3 轉輪葉片動應力對轉速的靈敏度曲線Fig.3 Sensitivity curve of blade dynamic stress to speed

4 結論

Sobol’法是一種基于方差的蒙特卡洛法，采用Sobol’法能夠有效地分析和計算轉輪葉片動態響應對工況參數的靈敏度。首先，建立轉輪葉片的動力學模型，推導轉輪葉片動態響應與各相關參數的內在關系式，并構建相關輸入參數的空間域；然后，將轉輪葉片的動態響應函數式表示為2k個 遞增項之和，并且使各項之間都是相互正交的；最后，基于轉輪葉片動態響應函數的總方差和偏方差表達式，得到轉輪葉片動態響應靈敏度計算公式。實例分析結果表明，轉輪葉片動態響應對轉速的靈敏度呈先增大后減小再增大的趨勢，在轉速∕額定轉速為0.6 以及1.2 時，轉輪葉片動態響應對轉速的靈敏度出現了兩個峰值，說明在這兩個轉速附近轉輪葉片動態響應對轉速的變化是較為敏感的。本文的研究為進一步深入分析轉輪葉片的振動特性與工況參數之間的內在關系提供一個可行的方法。