基于Rayleigh法的風(fēng)機塔筒預(yù)應(yīng)力基本頻率解析計算方法

雷振博， 劉 綱,2， 楊 微， 李 楊

(1.重慶大學(xué) 土木工程學(xué)院，重慶 400045; 2. 重慶大學(xué) 山地城鎮(zhèn)建設(shè)與新技術(shù)教育部重點實驗室，重慶 400045;3. 中國船舶重工集團(tuán)海裝風(fēng)電股份有限公司，重慶 401122; 4. 重慶大學(xué) 機械工程學(xué)院，重慶 400044)

風(fēng)電作為應(yīng)對生態(tài)環(huán)境污染和氣候變化的主力軍，逐漸得到世界各國的普遍共識和大力發(fā)展。統(tǒng)計表明，截至2019年我國風(fēng)電裝機容量已達(dá)2.21億kW，雄踞世界第一[1]。到2050年，我國風(fēng)電裝機容量更有望達(dá)到30億kW，屆時將滿足國內(nèi)30%以上的電力需求[2]。風(fēng)力發(fā)電技術(shù)將成為我國經(jīng)濟(jì)向綠色發(fā)展轉(zhuǎn)型、實現(xiàn)碳中和的重要手段。

設(shè)計過程中為避免出現(xiàn)共振且同時考慮經(jīng)濟(jì)安全等因素，塔筒第一階頻率(基頻)往往會遠(yuǎn)離風(fēng)輪旋轉(zhuǎn)頻率(1f)、過槳頻率(3f)以及主要的環(huán)境荷載頻率[3]，故對風(fēng)力機自振頻率的計算精度要求較高。而主流大型風(fēng)力機塔筒通常采用變截面錐形筒體形式[4]，截面剛度隨高度呈非線性變化，這給塔筒基頻的準(zhǔn)確計算帶來較大困難。

目前風(fēng)力機塔筒自振頻率計算主要有解析和數(shù)值計算兩類方法。在解析計算方面，Byrne[5]將塔筒結(jié)構(gòu)等效為懸臂梁，將機艙和風(fēng)輪等效為頂部集中質(zhì)量，基于結(jié)構(gòu)動力學(xué)原理，獲得塔筒等效單自由度體系基頻計算公式。牛文杰[6]考慮樁-彈性地基作用影響，采用柔度法對風(fēng)機結(jié)構(gòu)進(jìn)行多自由度動力分析，通過求解微分方程確定風(fēng)力機塔筒基頻。在數(shù)值計算方面，楊春寶等[7]考慮風(fēng)力機塔筒變截面特性，建立了近海風(fēng)力機塔筒自振頻率求解數(shù)值計算方法，并利用實際工程進(jìn)行有效性驗證。曾夢偉等[8]將風(fēng)輪與機艙等效為質(zhì)量點，建立機艙-塔筒-法蘭剛?cè)岫囿w耦合精細(xì)化有限元模型，采用數(shù)值模擬方法計算塔筒結(jié)構(gòu)自振頻率。田英鵬等[9]將風(fēng)力機法蘭和螺栓連接視為剛性連接，將風(fēng)輪、機艙及法蘭等效為質(zhì)量點，采用殼單元構(gòu)建塔筒結(jié)構(gòu)有限元模型進(jìn)行基頻計算，并通過試驗驗證了計算精度。

總體而言，采用精細(xì)化有限元模型能夠較為準(zhǔn)確獲得風(fēng)力機塔筒的基頻，但建模工作量極大，并不適合于塔筒前期設(shè)計的簡化、快速計算。而目前解析方法較為簡略，未考慮塔筒截面變剛度特性、機艙和風(fēng)輪載荷引起的塔筒預(yù)變形，計算精度不高，難以滿足工程實際需求。針對以上問題，根據(jù)風(fēng)機塔筒結(jié)構(gòu)構(gòu)造形式及力學(xué)特征，計入機艙及風(fēng)輪質(zhì)量對塔筒的預(yù)應(yīng)力作用，提出基于Rayleigh法的風(fēng)機塔筒預(yù)應(yīng)力基頻解析計算方法。同時，考慮塔筒變形形狀、頂部集中質(zhì)量偏心等多種情況，構(gòu)建塔筒基頻工程實用解析計算方法，為風(fēng)力機塔筒設(shè)計及共振控制提供方法支撐。

1 能量法基頻計算原理

1.1 風(fēng)力機塔筒模型簡化

錐臺型塔筒由數(shù)段鋼制錐筒通過法蘭盤連接而成，每段錐筒的直徑、壁厚均不一致，由底向上逐漸減小,如圖1(a)所示。大量研究、試驗和現(xiàn)場實測表明，塔筒主要沿順風(fēng)向和橫風(fēng)向振動，且往往以順風(fēng)向為主，故可將空間塔筒簡化為二維模型，取塔筒底面中心點為坐標(biāo)原點，沿豎直方向為z軸，沿順風(fēng)向為x軸，如圖1(b)所示。

圖1 風(fēng)機等效單自由度結(jié)構(gòu)Fig.1 Equivalent single degree of freedom flexible structure

由圖1可知，塔筒由n段鋼制塔段及n-1個法蘭盤組成，塔段從下至上的高度分別為l1，l2，…，ln，塔筒總高為L?？紤]塔筒沿高度方向具有變化的彎曲剛度EI(z)及質(zhì)量密度m(z)。法蘭盤質(zhì)量依下至上為m1，m2，…，mn-1，塔筒頂部機組及風(fēng)輪的質(zhì)量記為mn，其對塔筒的預(yù)應(yīng)力為集中質(zhì)量自質(zhì)量產(chǎn)生的軸向壓力N=mng。

1.2 能量法自振頻率計算原理

為方便采用能量法建立風(fēng)力機塔筒運動方程，可將塔筒近似等效為單自由度體系[10]。根據(jù)結(jié)構(gòu)動力學(xué)原理，假設(shè)塔筒某階振動的形狀函數(shù)為ψ(z)，則高度z處塔筒在t時刻的位移v(z,t)為

v(z,t)=ψ(z)Z(t)

(1)

式中：Z(t)為廣義坐標(biāo)，表示t時刻塔頂處x向的位移大小。

根據(jù)虛功原理可得塔筒廣義單自由度體系的運動方程為

(2)

其中，

(3)

式中：m*為廣義質(zhì)量；c*為廣義阻尼；k*為廣義彎曲剛度。

結(jié)合Rayleigh能量法和式(3)，可推導(dǎo)得到塔筒廣義單自由度體系的基頻為[11]

(4)

2 考慮預(yù)應(yīng)力的基頻解析方法

2.1 塔筒預(yù)變形

大量文獻(xiàn)表明，塔筒在頂部機艙、風(fēng)輪載荷作用下會產(chǎn)生壓應(yīng)力和預(yù)變形，塔筒中的壓應(yīng)力會導(dǎo)致塔筒剛度變化，進(jìn)而影響塔筒頻率[12]。因此，機艙、風(fēng)輪載荷引起壓預(yù)應(yīng)力效應(yīng)是影響風(fēng)機塔筒基頻的重要因素之一。

假定塔筒在機艙、風(fēng)輪載荷下產(chǎn)生的預(yù)變形為e(t)(見圖1(b)),根據(jù)幾何關(guān)系有

(5)

故考慮預(yù)應(yīng)力的塔筒基頻為

(6)

塔筒的廣義單自由度體系的運動方程為

(7)

(8)

式中，δe(t)為對預(yù)變形為e(t)的變分。

對比式(4)和式(6)可知，塔筒機艙、風(fēng)輪載荷作用下產(chǎn)生的壓應(yīng)力將降低塔筒基頻。

2.2 平均等效當(dāng)量慣性矩及壁厚

塔筒沿高度方向具有變化的壁厚及截面直徑，故簡化基頻計算時，采用等效壁厚t*來表示塔筒壁厚的影響[13]，其表達(dá)式為

(9)

式中：li為第i段塔筒高度；ti為第i段塔筒壁厚。

針對截面變直徑引起塔筒沿高度剛度變化問題，基于位移等效原理，采用平均等效當(dāng)量慣性矩I*近似計算[14]

(10)

式中：i1為系數(shù)，當(dāng)考慮塔筒頂部橫向集中荷載作用時i1=3，當(dāng)考慮塔筒橫向均布荷載作用時i1=4；Ii為通過積分獲得的第i段塔筒的慣性矩，其表達(dá)式為

(11)

式中：Di為第i段塔筒的底部中直徑；di為第i段塔筒的頂部中直徑。

3 多參數(shù)影響下基本頻率分析

3.1 考慮塔筒變形形狀影響

由Rayleigh法基本原理可知，基頻計算精度依賴于所取形狀函數(shù)與真實振型的相似度，對典型懸臂結(jié)構(gòu)，通常以集中荷載或均布荷載變形下的撓曲線作為形狀函數(shù)。風(fēng)力機塔筒結(jié)構(gòu)的力學(xué)模型，如圖2所示。

圖2 風(fēng)機變形影響因素Fig.2 Influence factors of fan deformation

3.1.1 塔頂集中荷載

如圖2(a)所示。根據(jù)材料力學(xué)可知，懸臂柱自由端作用水平集中荷載產(chǎn)生的撓曲線為[15]

(12)

將式(12)代入式(1)可得塔頂集中荷載下塔筒變形的形狀函數(shù)為

(13)

將式(13)代入式(6)可得集中荷載作用下，考慮預(yù)壓應(yīng)力影響的塔筒基頻表達(dá)式為

(14)

式中：ρ為材料密度；D為整個風(fēng)力機塔筒底部中直徑；d為頂部中直徑。

3.1.2 均布荷載

由材料力學(xué)可知，懸臂柱橫向作用水平均布荷載產(chǎn)生的撓曲線為

(15)

將式(15)代入式(1)可得均布荷載下塔筒變形的形狀函數(shù)為

(16)

將式(16)代入式(6)可得均布荷載作用下，考慮預(yù)壓應(yīng)力影響的塔筒基頻表達(dá)式為

(17)

3.2 考慮頂部質(zhì)量偏心影響

由于塔筒頂部的機艙、葉輪在x向有偏心，即圖2(b)中集中力N與塔筒的中心軸有偏心距Δ。塔頂質(zhì)量偏心將引起附加彎矩，則塔筒廣義幾何剛度應(yīng)修正為

(18)

將式(13)、式(18)代入式(6)，式(16)、式(18)代入式(6)可分別得到集中荷載作用下考慮偏心的塔筒基頻f3、均布荷載作用下考慮偏心的塔筒基頻f4為

(19)

4 數(shù)值算例

4.1 3.4 MW-100塔型

以兩個風(fēng)機塔筒為例，通過有限元數(shù)值計算結(jié)果驗證所提塔筒基頻計算公式的適用性。塔筒1為某3.4 MW風(fēng)機100 m高塔筒，該型風(fēng)機的機艙質(zhì)量為126.75 t，風(fēng)輪質(zhì)量為96.5 t；風(fēng)機塔筒整體呈圓臺狀，塔筒計算總高度為96.7 m；塔筒分節(jié)及截面詳細(xì)尺寸，如表1所示。塔筒材料采用Q345鋼，彈性模量取為2.06×105MPa，密度取為7 850 kg/m3。塔筒頂部風(fēng)輪沿x方向偏心-4.475 m，機艙沿x方向偏心0.953 m，沿y向偏心-0.038 m。

表1 塔筒幾何參數(shù)Tab.1 Tower geometry parameters

通過ANSYS軟件建立該塔筒的有限元殼模型，并計入塔頂機艙、風(fēng)輪質(zhì)量偏心的影響。塔筒選用SHELL63單元模擬，并考慮實際塔筒的漸變分段效果。頂部機艙和風(fēng)輪以及各法蘭盤等效為點質(zhì)量，采用MASS21單元模擬。塔筒底部與基礎(chǔ)固定連接，建立的有限元模型如圖3所示。

圖3 風(fēng)機塔筒有限元模型示意圖Fig.3 Finite element model of wind turbine tower

預(yù)應(yīng)力模態(tài)分析通過有預(yù)應(yīng)力的模態(tài)分析方法來實現(xiàn)，有預(yù)應(yīng)力結(jié)構(gòu)的固有頻率和模態(tài)計算，需要先進(jìn)行結(jié)構(gòu)的靜力分析獲得應(yīng)力剛度矩陣，將其代入頻率分析方程來得到考慮機艙和風(fēng)輪質(zhì)量預(yù)應(yīng)力的模態(tài)分析[16]。通過模態(tài)分析，獲取塔筒x向和y向的基頻如圖4所示。

從圖4可知，塔筒在x向和y向第一階頻率分別為0.222 Hz和0.225 Hz。根據(jù)式(14)、式(18)和式(19)，計算得到不同條件下的風(fēng)機塔筒基頻，如表2所示。

圖4 100 m風(fēng)機塔筒預(yù)應(yīng)力模態(tài)示意圖Fig.4 Pre-stress mode diagram of 100 m wind turbine tower

表2 100 m塔筒基頻計算結(jié)果對比Tab.2 Frequency calculation of 100 m wind turbine

4.2 3.4 MW-90塔型

塔筒2為某3.4 MW風(fēng)機90 m高塔筒。該型風(fēng)機的機艙質(zhì)量為126.75 t，風(fēng)輪質(zhì)量為101.72 t；風(fēng)機塔筒整體呈圓臺狀，塔筒計算總高度為86.7 m；塔筒頂部風(fēng)輪沿x方向偏心-4.475 m，機艙沿x方向偏心0.953 m,沿y向偏心-0.038 m。同4.1節(jié)相同建模方法，采用ANSYS有限元軟件計算得到塔筒在x向和y向基頻分別為0.262和0.265。根據(jù)式(14)、式(18)和式(19)，計算得到不同條件下的風(fēng)機塔筒基頻如表3所示。

表3 90 m塔筒頻率計算結(jié)果對比Tab.3 Frequency calculation of 90 m wind turbine

從表2、表3計算結(jié)果可知，無論x向還是y向，采用Rayleigh能量法所得基頻均高于精細(xì)化有限元模型計算值，這主要因為能量法計算所假設(shè)的變形公式式(13)、式(16)均無法與真實變形曲線完全吻合所致；相比而言，f2，f4較f1，f3更加接近于有限元數(shù)值解，因此，宜采用均布荷載作用下的形狀函數(shù)進(jìn)行解析法計算。同時，因x向偏心較大，考慮頂部質(zhì)量偏心對x向基頻計算精度有顯著影響，在集中荷載變形下，塔筒1和塔筒2的精度分別可提高2.8%和1.9%，均布荷載變形下分別可提高3.4%和2.6%。故當(dāng)頂部質(zhì)量偏心較大時，應(yīng)考慮偏心對解析法求解的影響。

5 試驗驗證

5.1 模型設(shè)計

以塔筒1為原型在實驗室搭建試驗塔筒，模型幾何相似比Sl為1/20。由于原型筒壁最薄處僅14 mm，無法完全遵守幾何相似關(guān)系進(jìn)行縮尺，故主要保證模型與原型動力特性相似，即保證構(gòu)件截面彎曲剛度相似與質(zhì)量分布相似[17-18]。通過材性試驗測得其彈性模量為1.525×105MPa，為保證加工精度可控，模型塔筒每節(jié)截面和壁厚不變，按原型塔筒相應(yīng)節(jié)底截面進(jìn)行相似設(shè)計，試驗塔筒模型幾何參數(shù)如表4所示。

表4 模型幾何參數(shù)Tab.4 Model geometric parameters

與原型塔筒對應(yīng)，模型塔筒總共分5段，段與段之間由法蘭連接，法蘭盤上設(shè)置配重籃，塔筒頂部設(shè)置配重箱考慮機艙和風(fēng)輪質(zhì)量。模型由Q345熱軋鋼板焊接而成組成，底部法蘭連接600 mm×600 mm×20 mm的矩形連接板，并通過18顆螺栓和鋼板連接。試驗?zāi)Ｐ晚敳颗渲叵溲豿方向偏心6 cm。筒體與法蘭之間采用對接焊縫連接，各法蘭間采用螺栓連接，如圖5所示。

5.2 測點及加載方案設(shè)計

在各配重籃和頂部配重箱的x和y向均布設(shè)P15H-2壓電加速度傳感器(見圖5)。

圖5 試驗?zāi)Ｐ虵ig.5 Experiment model

5.3 頻率測試

在塔頂施加錘擊激勵使其產(chǎn)生自由振動，測得塔頂加速度時程如圖6所示。通過傅里葉變換，得到塔筒兩個方向的功率譜如圖7所示。

圖6 塔頂加速度響應(yīng)時程曲線Fig.6 Acceleration response at the top of tower

圖7 加速度功率譜Fig.7 Acceleration power spectrum

從圖7可看出，模型對應(yīng)x向和y向基頻分別為1.18 Hz和1.22 Hz。同時，試驗中僅采用頂部激勵時，圖6的加速度響應(yīng)出現(xiàn)典型的拍頻振動特征及正交耦合現(xiàn)象，表明試驗?zāi)Ｐ偷恼駝臃细呷峤Y(jié)構(gòu)的響應(yīng)特征[19-21]。加速度功率譜僅有一個峰值，表明塔筒以第一階振動為主，這與大多數(shù)文獻(xiàn)及試驗結(jié)果相吻合[22-24]。

5.4 有限元數(shù)值模型

與4.1節(jié)相同建模方法相同，彈性模量采用材性試驗測得的1.525×105MPa。結(jié)合表4模型幾何參數(shù)，采用ANSYS有限元軟件計算試驗縮尺模型的x向和y向基頻分別為1.203 Hz和1.212 Hz。x向和y向的基本頻率與試驗測試的基頻誤差均在3%以內(nèi)，與試驗測試值吻合較好。

5.5 解析法對比驗證

結(jié)合傳統(tǒng)風(fēng)力機塔筒頻率解析方法及式(18)，計算可得試驗?zāi)Ｐ徒馕龇ǖ睦碚撝蹬c動力測試基頻，如表5所示。

表5 理論值與實測值對比Tab.5 Comparison between calculation and experiment

從表5可知，試驗塔筒在x向?qū)崪y基頻低于y向基頻，這主要是因為該塔筒在x向存在偏心，將減少塔筒在該向剛度，故實測頻率值偏小，故在計算塔筒基頻時有必要計入偏心影響。就傳統(tǒng)解析法和本文所提方法而言，后者較前者更加精確，即x向計算精度提高4.52%，y向計算精度提高4.50%。這主要是因為本文方法綜合考慮了機艙、風(fēng)輪質(zhì)量和偏心，以及塔筒不同變形形狀等情況的影響。

6 結(jié) 論

針對風(fēng)力機塔筒基頻解析計算精度難以滿足工程要求的問題，提出一種基于Rayleigh法的風(fēng)力機塔筒基頻快速估算方法?？紤]機艙、風(fēng)輪質(zhì)量對塔筒壓力作用，以及塔筒不同變形形狀函數(shù)、頂部集中質(zhì)量偏心等系列實際情況，基于理論推導(dǎo)給出了解析計算公式，并通過數(shù)值算例與縮尺模型試驗驗證了所提方法的適用性及有效性。所得結(jié)論如下：

(1) 有限元模型與試驗結(jié)果均表明采用均布荷載作用下的變形形狀曲線將得到更高的解析計算精度，兩個數(shù)值算例中，f2和f4與有限元模型計算值的最大相對誤差僅為3.47%。因此，宜采用均布荷載作用下的撓曲線作為Rayleigh法的形狀函數(shù)。

(2) 頂部質(zhì)量偏心對塔筒基頻計算精度有顯著影響，特別是在偏心距較大時。采用集中荷載作用下的變形曲線時，考慮偏心情況下塔筒1和塔筒2的精度分別可提高2.8%和1.9%；采用均布荷載作用下變形曲線時，考慮偏心情況下塔筒1和塔筒2的計算精度分別提高3.4%和2.6%。因此，在計算風(fēng)機塔筒基頻時，應(yīng)考慮機艙、風(fēng)輪質(zhì)量偏心的影響。

(3) 縮尺模型試驗表明，本文方法考慮了機艙、風(fēng)輪質(zhì)量及偏心，以及塔筒不同變形形狀等多種情況的綜合影響，x向基頻計算精度可提高4.52%，y向基頻計算精度可提高4.50%，較傳統(tǒng)解析方法更加精確。