風(fēng)力機(jī)塔架－葉片耦合模型風(fēng)致響應(yīng)時(shí)域分析＊

柯世堂，曹九發(fā)，王 瓏，王同光

（南京航空航天大學(xué) 江蘇省風(fēng)力機(jī)設(shè)計(jì)高技術(shù)研究重點(diǎn)實(shí)驗(yàn)室，江蘇 南京 210016）

風(fēng)荷載是風(fēng)力機(jī)系統(tǒng)設(shè)計(jì)最重要荷載之一[1].隨著兆瓦級(jí)大功率風(fēng)力機(jī)的普遍應(yīng)用，葉片的旋翼直徑達(dá)到百米量級(jí)，從而滿足捕捉風(fēng)能的要求，但同時(shí)也遭受了強(qiáng)大風(fēng)推力，并通過(guò)機(jī)艙將風(fēng)荷載傳遞給塔架，使得葉片和塔架之間的耦合效應(yīng)愈加明顯.以往針對(duì)風(fēng)力機(jī)系統(tǒng)的抗風(fēng)研究，關(guān)注風(fēng)力機(jī)空氣動(dòng)力學(xué)性能的研究者大多忽略塔架的影響，對(duì)葉片的流場(chǎng)特性[2－3]和氣彈響應(yīng)[4－5]進(jìn)行了深入詳細(xì)的研究；而關(guān)注風(fēng)力機(jī)結(jié)構(gòu)受力性能的研究者大多以塔架為研究對(duì)象，把葉片遭受的風(fēng)載等效為擬靜力施加到塔頂進(jìn)行風(fēng)振計(jì)算，或是建模時(shí)考慮塔架－葉片的耦合模型，而風(fēng)場(chǎng)模擬時(shí)沒(méi)有完全考慮葉片平穩(wěn)風(fēng)修正、旋轉(zhuǎn)效應(yīng)和氣彈效應(yīng)等影響[6－8].

目前，對(duì)風(fēng)力機(jī)全機(jī)結(jié)構(gòu)的抗風(fēng)設(shè)計(jì)，基本上采用等效風(fēng)荷載的擬靜力分析方法，少量文獻(xiàn)進(jìn)行了頻域內(nèi)的風(fēng)振分析，考慮到塔架－葉片耦合模型脈動(dòng)風(fēng)場(chǎng)模擬的復(fù)雜性，在時(shí)域內(nèi)對(duì)風(fēng)力機(jī)塔架－葉片耦合模型進(jìn)行風(fēng)致響應(yīng)計(jì)算則更為復(fù)雜.然而時(shí)域分析可以克服頻域分析中基于線性化假設(shè)的不足，同時(shí)時(shí)域計(jì)算結(jié)果的直觀性有利于工程師了解風(fēng)力機(jī)系統(tǒng)的受力性能.因此，基于風(fēng)力機(jī)塔架－葉片耦合模型的脈動(dòng)風(fēng)荷載模擬結(jié)果，基于塔架－葉片耦合有限元模型進(jìn)行風(fēng)振響應(yīng)時(shí)域分析，對(duì)風(fēng)力機(jī)系統(tǒng)抗風(fēng)的精細(xì)化設(shè)計(jì)有重要的指導(dǎo)作用和實(shí)用價(jià)值.

鑒于此，本文基于時(shí)域方法分析風(fēng)力機(jī)塔架－葉片耦合模型在順風(fēng)向風(fēng)荷載動(dòng)力作用下的響應(yīng)特性，建立了風(fēng)力機(jī)塔架－葉片耦合的計(jì)算模型，運(yùn)用諧波合成法和葉素－動(dòng)量理論進(jìn)行數(shù)值模擬獲得了考慮平穩(wěn)風(fēng)修正、葉片旋轉(zhuǎn)效應(yīng)和相干性的塔架－葉片耦合模型氣動(dòng)載荷，在此基礎(chǔ)上利用有限元方法進(jìn)行了風(fēng)力機(jī)系統(tǒng)風(fēng)致響應(yīng)時(shí)域計(jì)算，基于計(jì)算結(jié)果，分析了風(fēng)力機(jī)塔架－葉片耦合模型的風(fēng)致響應(yīng)特性.

1 風(fēng)力機(jī)塔架－葉片耦合模型風(fēng)場(chǎng)模擬

由于風(fēng)力機(jī)塔架－葉片耦合系統(tǒng)結(jié)構(gòu)本身的復(fù)雜性及隨機(jī)脈動(dòng)風(fēng)作用的不確定性，風(fēng)場(chǎng)模擬時(shí)需作如下簡(jiǎn)化：①只考慮順風(fēng)向風(fēng)荷載數(shù)值模擬和作用；②脈動(dòng)風(fēng)被假定為零均值的平穩(wěn)高斯隨機(jī)過(guò)程；③不考慮葉片和氣流之間的氣動(dòng)彈性作用.

1.1 平穩(wěn)風(fēng)修正

風(fēng)力機(jī)的平穩(wěn)風(fēng)速由于受到風(fēng)切變、塔影效應(yīng)和上游風(fēng)機(jī)尾流的影響，風(fēng)場(chǎng)模擬時(shí)必須要對(duì)平穩(wěn)風(fēng)模型進(jìn)行修正[9].其中風(fēng)剪切主要采用指數(shù)模型，塔影效應(yīng)主要采用適用于葉片在塔架上風(fēng)向運(yùn)行的潛流模型，上游尾流影響主要采用包含由于尾流引起的附加湍流的渦流粘度模型.

風(fēng)剪切的影響是指平穩(wěn)風(fēng)速隨著高度而變化，其指數(shù)修正模型表達(dá)公式如下：

式中：V（h）為指高度h處的風(fēng)速；V（h0）為輪轂的參考高度h0處的參考風(fēng)速；h0為輪轂的位置；當(dāng)不考慮風(fēng)剪切的影響時(shí)，可以將α的值設(shè)為0，取值一般為0.1～0.25.

塔影是由于風(fēng)力機(jī)塔架的存在影響了風(fēng)場(chǎng)平穩(wěn)風(fēng)速，其影響修正主要有3種模型：葉片在塔架上風(fēng)向運(yùn)行的潛流模型、葉片在塔架下風(fēng)向運(yùn)行的經(jīng)驗(yàn)?zāi)Ｐ秃蛢烧叩慕M合模型.現(xiàn)有風(fēng)力機(jī)基本都是上風(fēng)向風(fēng)機(jī)設(shè)計(jì)，用以避免葉片周期地通過(guò)塔架尾跡產(chǎn)生的附加噪聲和激振力，因此采用潛流模型修正：

式中：DT為開(kāi)始考慮塔影影響的高度處的塔架直徑；F為塔架直徑修正因子；z為計(jì)算點(diǎn)到塔架中心的縱向距離；x為風(fēng)矢量經(jīng)過(guò)時(shí)距離塔架中心橫向距離.

由于風(fēng)力機(jī)很少會(huì)存在單個(gè)運(yùn)營(yíng)狀態(tài)，基本都是以群體形式建立.因此，在模擬風(fēng)力機(jī)模型的風(fēng)場(chǎng)時(shí)，當(dāng)由于單個(gè)風(fēng)力機(jī)系統(tǒng)轉(zhuǎn)子部分或者全部處于上游風(fēng)力機(jī)的尾流中時(shí)，就必須考慮上游風(fēng)力機(jī)尾流對(duì)風(fēng)力機(jī)風(fēng)場(chǎng)的影響.本文采用渦流粘度模型來(lái)考慮尾流引起的附加湍流影響，修正公式如下：

式中：r為局部速度到尾流中心的距離；w為尾流寬度；V0為未受擾動(dòng)平均風(fēng)速；Δ為速度逆差.

1.2 脈動(dòng)風(fēng)模擬方法

已有研究表明，風(fēng)力機(jī)葉片的風(fēng)場(chǎng)模擬可采用葉素－動(dòng)量理論[10]，相較渦尾跡方法和 CFD方法，葉素－動(dòng)量理論可以實(shí)現(xiàn)葉片旋轉(zhuǎn)風(fēng)場(chǎng)的模擬，且計(jì)算量較小、耗時(shí)較短，滿足風(fēng)力機(jī)日?？焖儆?jì)算的要求.考慮到葉素－動(dòng)量理論是針對(duì)葉片在某一特定來(lái)流風(fēng)速和工況下對(duì)應(yīng)參考風(fēng)速的求解，本文通過(guò)諧波合成法[11]模擬獲得考慮葉片和塔架相干性的來(lái)流風(fēng)速時(shí)程，再對(duì)每個(gè)來(lái)流風(fēng)速樣本采用改進(jìn)的葉素－動(dòng)量理論進(jìn)行計(jì)算得到該來(lái)流風(fēng)速下葉片的參考風(fēng)速，從而考慮葉片的旋轉(zhuǎn)效應(yīng)和相干性，如此迭代循環(huán)獲得葉片的參考風(fēng)速時(shí)程.風(fēng)力機(jī)塔架－葉片耦合模型脈動(dòng)風(fēng)模擬流程示意圖見(jiàn)圖1.

圖1 風(fēng)力機(jī)塔架－葉片耦合模型風(fēng)場(chǎng)模擬流程示意圖Fig.1 The flow sketch figure of wind field simulation for wind turbine tower－blade coupled model

風(fēng)譜模型采用改進(jìn)的Von karman模型，因?yàn)橄啾菵avenport模型和Kaimal模型，改進(jìn)的Von karman風(fēng)譜模型能夠更加準(zhǔn)確地產(chǎn)生隨機(jī)風(fēng)場(chǎng)樣本，該模型校正了基本模型在高度150m以下的缺陷，更加符合風(fēng)力機(jī)自身風(fēng)場(chǎng)特性.圖2為3種風(fēng)譜模型曲線，在低頻區(qū)3種風(fēng)譜數(shù)值接近，在峰值附近區(qū)域，Kaimal和Von karman頻譜分別要低于和高于改進(jìn)Von karman頻譜；而在高頻區(qū)域，Kaimal和Von karman頻譜幾乎重合，但是改進(jìn)Von karman頻譜會(huì)略高于另外兩個(gè)頻譜幅值.

圖2 3種標(biāo)準(zhǔn)化風(fēng)譜曲線示意圖Fig.2 The different longitudinal normal powerspectrum

根據(jù)風(fēng)力機(jī)風(fēng)場(chǎng)的自身特性，選擇改進(jìn)的Von karman風(fēng)譜模型：

式中：Suu為隨機(jī)風(fēng)速變化自頻譜；f為頻率；σu為隨機(jī)風(fēng)的標(biāo)準(zhǔn)差；~nu為無(wú)量綱頻率參數(shù)＝fLux／U10，U10為10m高基本風(fēng)速，Lux為隨機(jī)風(fēng)縱向分量長(zhǎng)度尺度.

采用Davenport相關(guān)系數(shù)經(jīng)驗(yàn)公式來(lái)考慮葉片和塔架之間的相關(guān)性：

式中：Cx，Cy和Cz分別為葉片和塔架上任意2點(diǎn)橫向、順風(fēng)向和垂直向的衰減系數(shù)；ω為脈動(dòng)風(fēng)頻率；v（H）為H處平均風(fēng)速.僅考慮順風(fēng)向相關(guān)性時(shí)，其余兩個(gè)方向衰減系數(shù)均取0.

定義風(fēng)機(jī)上n個(gè)風(fēng)速模擬節(jié)點(diǎn)，假定均為零均值的平穩(wěn)高斯過(guò)程，其風(fēng)譜密度函數(shù)矩陣為：

式中：sii（ω）為節(jié)點(diǎn)脈動(dòng)風(fēng)自功率譜，采用式（2）中的風(fēng)譜模型計(jì)算；sij（ω）為互功率譜，其表達(dá)式需要用到塔架和葉片以及塔架上任意兩點(diǎn)之間的相干性，計(jì)算公式為：

其中葉片需考慮旋轉(zhuǎn)平面內(nèi)的各點(diǎn)間的相干性，以及葉片和塔架之間的相干性影響.再將S（ω）進(jìn)行Cholesky分解，此時(shí)風(fēng)力發(fā)電塔架上的任何一個(gè)節(jié)點(diǎn)脈動(dòng)風(fēng)速時(shí)程可以由其功率譜決定，根據(jù)Shinozuka理論，模擬的風(fēng)速時(shí)程可以表達(dá)為：

式中：風(fēng)譜在頻率范圍內(nèi)劃分成N個(gè)相同部分，Δω＝ω／N為頻率增量；｜Hjm（ωl）｜為基于 Davenport來(lái)流風(fēng)譜矩陣進(jìn)行Cholesky分解獲得的下三角矩陣的模；θml為介于0和2π之間均勻分布的隨機(jī)數(shù)，可采用Matlab的隨機(jī)數(shù)生成函數(shù)，建議每次生成隨機(jī)數(shù)后應(yīng)恢復(fù)初始狀態(tài)；ωl＝lΔω為頻域的遞增變量；ψjm（ωl）為兩個(gè)不同作用點(diǎn)之間的相位角，它是由Hjm（ωl）的虛部和實(shí)部的比值確定.

采用修正的BEM 理論[12]，引入葉根損失和葉尖損失，在軸向誘導(dǎo)因子較大時(shí)使用Ct的經(jīng)驗(yàn)?zāi)Ｐ?，并加入?dòng)態(tài)入流和動(dòng)態(tài)失速模型.使用該方法，可以計(jì)算風(fēng)力機(jī)在不同風(fēng)速、轉(zhuǎn)速、槳距角及偏航角情況下的動(dòng)態(tài)載荷.進(jìn)而獲取作用在葉片上的脈動(dòng)風(fēng)速時(shí)程.具體方法如下：

根據(jù)BEM理論，葉片上的相對(duì)風(fēng)速Vrel采用下式計(jì)算：

式中：vox，voy分別為沿順風(fēng)向和橫風(fēng)向的來(lái)流脈動(dòng)風(fēng)速，采用式（9）諧波疊加法計(jì)算；vrot為葉片旋轉(zhuǎn)導(dǎo)致的線速度；W為誘導(dǎo)速度；vbx和vby為葉片振動(dòng)速度.

誘導(dǎo)速度W可由下式表示：

式中：B為葉片數(shù)；L為指升力；φ為入流角；ρ為空氣密度；r為葉片截面的展向位置；n為推力方向的單位向量；F為普朗特葉尖損失因子；fg為Glauert修正.同時(shí)，本文還采用了動(dòng)態(tài)入流模型和動(dòng)態(tài)失速模型，修正葉片運(yùn)轉(zhuǎn)的非定常效應(yīng).

根據(jù)下式計(jì)算葉片攻角α：

式中：β為槳矩角；θtwist為葉片剖面幾何扭角，計(jì)算公式為：

1.3 數(shù)值模擬結(jié)果

本文以某5MW三槳葉變槳距特大型風(fēng)力機(jī)系統(tǒng)為例進(jìn)行風(fēng)場(chǎng)數(shù)值模擬和風(fēng)致響應(yīng)時(shí)域分析.塔高124m，底徑4.8m，頂徑2.6m，塔體通長(zhǎng)為變厚度結(jié)構(gòu)，底壁厚150mm，頂壁厚60mm，通長(zhǎng)厚度由底部至頂部呈線性減小趨勢(shì).機(jī)艙長(zhǎng)12m，寬4.6m，高4.2m，總質(zhì)量140.2×103kg.各槳葉間成120°夾角，沿周向平均分布，風(fēng)輪直徑為120m，寬度2.4m，厚度0.38m，長(zhǎng)度60m，偏航角為0°，額定轉(zhuǎn)速為17r／min.

基于上述方法，采用Matlab語(yǔ)言編制相應(yīng)的數(shù)值模擬程序MBEM－HSM.在計(jì)算過(guò)程中脈動(dòng)風(fēng)上限頻率取為2π，脈動(dòng)風(fēng)頻率分割點(diǎn)數(shù)取為2 048，頻率增量Δω＝0.003 07Hz.圖3和4給出了10m高平均風(fēng)速為24m／s模擬得到的葉片和塔架中部的風(fēng)速時(shí)程曲線和對(duì)應(yīng)功率譜曲線，為對(duì)比顯示風(fēng)力機(jī)塔架－葉片耦合模型風(fēng)場(chǎng)的特殊性，風(fēng)譜坐標(biāo)系均為對(duì)數(shù)坐標(biāo).

圖3 葉片參考風(fēng)速時(shí)程模擬結(jié)果Fig.3 Simulating result of fluctuating wind velocity of blades

圖4 塔架中部風(fēng)速時(shí)程模擬結(jié)果Fig.4 Simulating result of fluctuating wind velocity of tower middle part

對(duì)比發(fā)現(xiàn)，考慮塔架－葉片耦合模型葉片脈動(dòng)風(fēng)速功率譜曲線在高頻處存在明顯的能量浮動(dòng)，應(yīng)該是由于考慮改進(jìn)的Von karman風(fēng)譜模型、三維旋轉(zhuǎn)和干擾效應(yīng)而產(chǎn)生的現(xiàn)象，在風(fēng)力機(jī)系統(tǒng)風(fēng)場(chǎng)模擬中應(yīng)引起重視；塔架中部的風(fēng)速功率譜和Davenport風(fēng)譜吻合較好，這是由于隨著高度的降低，塔架受到葉片的干擾效應(yīng)減小，在高頻出現(xiàn)的能量浮動(dòng)和數(shù)值逐漸減小.根據(jù)本文模擬過(guò)程和對(duì)比分析可以認(rèn)為，采用本文的MBEM－HSM方法可以很好地模擬風(fēng)力機(jī)塔架－葉片耦合模型風(fēng)場(chǎng).

2 風(fēng)力機(jī)塔架－葉片耦合模型動(dòng)力特性

基于ANSYS軟件平臺(tái)，建立了風(fēng)力機(jī)塔架－葉片耦合模型.其中葉片和塔體采用SHELL91單元，機(jī)艙及其內(nèi)部結(jié)構(gòu)可作為整體采用梁?jiǎn)卧狟EAM189模擬.通過(guò)多點(diǎn)約束單元耦合命令將各部分連接在一起，形成整體的風(fēng)力機(jī)塔架－葉片耦合模型.依據(jù)效率和精度均衡的原則，模型一共劃分了3 812個(gè)單元.動(dòng)力特性分析時(shí)把葉片旋轉(zhuǎn)產(chǎn)生的離心力作為預(yù)應(yīng)力預(yù)先均勻施加在葉片上，后續(xù)計(jì)算的頻率和模態(tài)信息均考慮葉片轉(zhuǎn)動(dòng)帶來(lái)的離心力作用.圖5給出了考慮／不考慮葉片離心力作用2種工況下系統(tǒng)前200階自振頻率的分布情況.

圖5 風(fēng)力機(jī)塔架－葉片耦合系統(tǒng)固有頻率分布圖Fig.5 Scattergram of natural frequence for wind turbine tower－blade coupled system

從圖5可以看出：1）風(fēng)力機(jī)塔架－葉片耦合模型的基頻很小，僅為0.24Hz，第50階自振頻率為9.2 Hz，各模態(tài)之間的間隔很小，其自振頻率數(shù)值遠(yuǎn)遠(yuǎn)小于相同高度和直徑的高聳結(jié)構(gòu)，說(shuō)明葉片和塔架的耦合效應(yīng)明顯降低了風(fēng)力機(jī)系統(tǒng)的基頻，使得系統(tǒng)結(jié)構(gòu)可能會(huì)對(duì)風(fēng)荷載更為敏感.2）當(dāng)考慮葉片轉(zhuǎn)動(dòng)引起的離心力作用時(shí)，耦合系統(tǒng)的基頻（0.27Hz）要略大于不考慮離心力作用下的系統(tǒng)基頻，并且隨著模態(tài)數(shù)目的增加，離心力效應(yīng)帶來(lái)的頻率影響越來(lái)越大.這是由于風(fēng)力機(jī)塔架－葉片耦合系統(tǒng)的高階頻率主要以葉片的左右擺動(dòng)和前后振動(dòng)為主，此時(shí)葉片旋轉(zhuǎn)引起的離心力效應(yīng)更為明顯，相應(yīng)增強(qiáng)了風(fēng)力機(jī)系統(tǒng)的剛度.本文的后續(xù)計(jì)算均采用考慮葉片離心力作用的更加真實(shí)的模態(tài)參數(shù).

從圖6給出的風(fēng)力機(jī)塔架－葉片系統(tǒng)模態(tài)振型圖中可以看出，塔架的彎曲變形與葉片的揮舞／擺振相互耦合.通過(guò)多階模型分析表明，系統(tǒng)振型主要以葉片的前后揮舞和左右擺動(dòng)為主，高階模態(tài)下塔架本身也會(huì)出現(xiàn)結(jié)構(gòu)變形和失穩(wěn)形態(tài)，并且和葉片變形耦合作用.動(dòng)力特性分析表明，風(fēng)力機(jī)系統(tǒng)動(dòng)力響應(yīng)計(jì)算必須要考慮葉片和塔架耦合作用.

圖6 風(fēng)力發(fā)電塔－輪系統(tǒng)的典型模態(tài)振型示意圖Fig.6 Typical modes of vibration on wind rotors systems

3 風(fēng)致響應(yīng)時(shí)域分析方法

根據(jù)達(dá)朗貝爾原理，風(fēng)力機(jī)系統(tǒng)結(jié)構(gòu)風(fēng)致動(dòng)力響應(yīng)控制方程為：

式中：p（t）為風(fēng)荷載向量；M，C和K分別為結(jié)構(gòu)質(zhì)量、阻尼和剛度矩陣;（t），（t）和y（t）分別為節(jié)點(diǎn)加速度向量、速度向量和位移向量.

通過(guò)將節(jié)點(diǎn)位移向量從物理坐標(biāo)轉(zhuǎn)換到廣義模態(tài)位移后，實(shí)現(xiàn)方程解耦，再引入模態(tài)阻尼，則各廣義模態(tài)對(duì)應(yīng)的運(yùn)動(dòng)方程為：

式中：qi（t）為第i階廣義模態(tài)位移響應(yīng)向量；mi，ci，ki和pi（t）分別為第i階模態(tài)的模態(tài)質(zhì)量、模態(tài)阻尼、模態(tài)剛度和模態(tài)力.這樣單自由度的運(yùn)動(dòng)方程可根據(jù)Duhamel積分原理，初始條件為零，數(shù)值解為：

式中：ωni＝（ki／mi）1／2為結(jié)構(gòu)的固有模態(tài)頻率；ζi＝bi／（2miωni）為模態(tài)阻尼比；ωdi＝ωni（1－ζ2）1／2為結(jié)構(gòu)阻尼振動(dòng)頻率；Δτ為積分時(shí)間步長(zhǎng).本文選用采用Newmark－β逐步積分法和Newton－Raphson迭代理論，其逐步積分方法的基本假設(shè)為：

再將各模態(tài)下的廣義位移轉(zhuǎn)換為物理位移并進(jìn)行疊加，可得到各節(jié)點(diǎn)風(fēng)致響應(yīng)位移為：

式中：Φ為系統(tǒng)振型矩陣.

風(fēng)力機(jī)塔架－葉片耦合模型風(fēng)致響應(yīng)時(shí)域計(jì)算基于有限元軟件ANSYS平臺(tái)，將第2節(jié)數(shù)值模擬得到的風(fēng)荷載時(shí)程作為外部激勵(lì)作用于風(fēng)力機(jī)塔架－葉片耦合模型上，采用Newmark－β逐步積分法和Newton－Raphson迭代理論，其中各模態(tài)阻尼均為0.02，積分時(shí)間步長(zhǎng)取為0.05s，響應(yīng)輸出時(shí)間步長(zhǎng)取為0.025s，截取模態(tài)取系統(tǒng)的前50階.

4 風(fēng)力機(jī)塔架－葉片耦合模型風(fēng)致響應(yīng)

當(dāng)平均風(fēng)速過(guò)大時(shí)，風(fēng)力發(fā)電機(jī)組會(huì)自動(dòng)停止運(yùn)行，這個(gè)臨界風(fēng)速稱之為切出風(fēng)速，對(duì)于兆瓦級(jí)風(fēng)力機(jī)通常設(shè)定為25m／s左右.首先，進(jìn)行切出風(fēng)速下風(fēng)力機(jī)塔架－葉片耦合模型風(fēng)致響應(yīng)時(shí)域計(jì)算.限于篇幅，本文只給出了風(fēng)力機(jī)耦合模型葉片尖部、機(jī)艙和塔架頂部位移響應(yīng)時(shí)程曲線，如圖7所示.

對(duì)比發(fā)現(xiàn)，基于風(fēng)力機(jī)塔架－葉片耦合模型時(shí)域計(jì)算得到的葉片尖部位移響應(yīng)以低頻共振效應(yīng)為主，其脈動(dòng)覆蓋范圍和數(shù)值最大，這是由于風(fēng)力機(jī)遭受的風(fēng)荷載卓越頻率與葉片的擺振／舞動(dòng)頻率較為接近，從而引發(fā)了葉片的共振.機(jī)艙由于剛度較大，風(fēng)荷載作用主要導(dǎo)致靜力或準(zhǔn)靜力響應(yīng)，因此位移響應(yīng)時(shí)程曲線浮動(dòng)范圍較小，均值約在0.295m.塔架頂部的位移響應(yīng)時(shí)程明顯小于葉片響應(yīng)數(shù)值，且高頻振型貢獻(xiàn)更加明顯，對(duì)應(yīng)的振型大多為葉片振動(dòng)和塔架變形的耦合作用.分析表明，考慮了葉片的風(fēng)場(chǎng)影響和結(jié)構(gòu)模態(tài)耦合作用后，塔架的風(fēng)致響應(yīng)特性明顯呈現(xiàn)多振型響應(yīng)和多荷載形態(tài).

圖7 風(fēng)力機(jī)塔架－葉片耦合模型典型部位位移響應(yīng)時(shí)程曲線Fig.7 The displacement time history on typical nodes of wind turbine tower－blade coupled model

為研究葉片對(duì)風(fēng)力機(jī)塔架風(fēng)致響應(yīng)特性的影響，再對(duì)不同風(fēng)速下考慮和不考慮葉片影響塔架風(fēng)致響應(yīng)進(jìn)行分析.塔頂順風(fēng)向位移和塔底彎矩脈動(dòng)響應(yīng)均方根隨風(fēng)速的變化曲線如圖8和圖9所示，由圖可見(jiàn)，在低風(fēng)速（≤10m／s）下是否考慮葉片影響對(duì)風(fēng)力機(jī)塔架風(fēng)致響應(yīng)影響不大，而隨著風(fēng)速的逐漸增大，考慮葉片影響計(jì)算得到的塔架風(fēng)致響應(yīng)要遠(yuǎn)大于不考慮葉片作用時(shí)的計(jì)算結(jié)果，其中考慮葉片影響時(shí)塔架頂部位移較未考慮的要大將近2倍，對(duì)于塔架底部彎矩來(lái)說(shuō)要大1.5倍.

圖8 不同風(fēng)速下風(fēng)力機(jī)塔架頂部順風(fēng)向位移響應(yīng)均方差Fig.8 RMS of along－wind displacement responses of tower top under different wind speed

圖9 不同風(fēng)速下風(fēng)力機(jī)塔底彎矩響應(yīng)均方差Fig.9 RMS of bending moment responses of tower bottom under different wind speed

5 結(jié)束語(yǔ)

基于風(fēng)力機(jī)塔架－葉片耦合模型，實(shí)現(xiàn)了給定風(fēng)環(huán)境的系統(tǒng)風(fēng)場(chǎng)模擬和考慮離心力效應(yīng)的動(dòng)力特性分析，并在此基礎(chǔ)上進(jìn)行了系統(tǒng)風(fēng)致響應(yīng)的時(shí)域分析.

數(shù)值算例分析表明，風(fēng)力機(jī)葉片風(fēng)場(chǎng)存在高頻能量，塔架風(fēng)譜更接近自然風(fēng)譜；塔架－葉片耦合模型以塔架變形伴隨著葉片振動(dòng)耦合振型為主，且葉片旋轉(zhuǎn)引起的離心力作用不能忽略；葉片風(fēng)致響應(yīng)以低階模態(tài)的共振效應(yīng)為主，塔架的風(fēng)致響應(yīng)中高階模態(tài)的貢獻(xiàn)逐漸顯著，忽略葉片影響帶來(lái)的誤差非常顯著.綜上所述，風(fēng)力機(jī)的抗風(fēng)設(shè)計(jì)應(yīng)當(dāng)考慮塔架－葉片的耦合效應(yīng)，本文仿真算法可為風(fēng)力機(jī)系統(tǒng)風(fēng)場(chǎng)模擬和風(fēng)致響應(yīng)時(shí)域分析提供參考.

[1] BUTTERFIELD S，MUSIAL W，SCOTT G.Definition of a 5－MW reference wind turbine for offshore system development[M].Colorado：National Renewable Energy Laboratory，2009：129－142.

[2] BAZILEVS Y，HSU M C，KIENDL J，etal.A computational procedure for prebending of wind turbine blades[J].International Journal for Numerical Methods in Engineering，2012，89（3）：323－336.

[3] HOOGEDOORN E，JACOBS G B，BEYENE A.Aero－elastic behavior of a flexible blade for wind turbine application：a 2D computational study[J].Energy，2010，35（2）：778－785.

[4] 金志昊，范宣華，蘇先樾，等.風(fēng)力機(jī)葉片順風(fēng)向風(fēng)致振動(dòng)研究[J].南京航空航天大學(xué)學(xué)報(bào)，2011，43（5）：677－681.JIN Zhi－h(huán)ao，F(xiàn)AN Xuan－h(huán)ua，SU Xian－yu，etal.Research on blade wind－induced vibration in wind direction [J].Journal of Nanjing University of Aeronautics ＆ Astronautics，2011，43（5）：677－681.（In Chinese）

[5] RONOLD K O，LARSEN G C.Optimization of a design code for wind turbine rotor blades in fatigue [J].Engineering Structure，2001，23：993－1002.

[6] 李德源，劉勝祥.風(fēng)波聯(lián)合作用下的風(fēng)力機(jī)塔架疲勞特性分析[J].太陽(yáng)能學(xué)報(bào)，2009，30（10）：1250－1256.LI De－yuan，LIU Sheng－xiang.Fatigue characteristics analysis of wind turbine tower under wind－wave combined effect[J].Acta Energiae Solaris Sinica，2009，30（10）：1250－1256.（In Chinese）

[7] NAESS A，GAIDAI O，HAVER S.Efficient estimation of extreme response of drag－dominated offshore structures by monte carlo simulation [J].Ocean Engineering，2007，34：2188－2197.

[8] 柯世堂，王同光，趙林，等.風(fēng)力機(jī)風(fēng)振背景、共振響應(yīng)特性及耦合項(xiàng)分析[J].中國(guó)電機(jī)工程學(xué)報(bào)，2013，33（26）：101－108.KE Shi－tang，WANG Tong－guang，ZHAO Lin，etal.Background，resonant components and coupled effect of wind－induced responses on wind turbine systems[J].Proceeding of the CSEE，2013，33（26）：101－108.（In Chinese）

[9] 廖明夫.風(fēng)力發(fā)電技術(shù)[M].西安：西北工業(yè)大學(xué)出版社，2008：60－72.LIAO Ming－fu.Wind turbine technologies[M].Xi’an：The Northwest Industry University Press，2008：60－72.（In Chinese）

[10] 伍艷，謝華，王同光.風(fēng)力機(jī)葉片的非定常氣動(dòng)特性計(jì)算方法的改進(jìn)[J].工程力學(xué)，2008，25（10）：54－60.WU Yan，XIE Hua，WANG Tong－guang.Modification of calculating unsteady aerodynamic characteristics of wind turbine blades[J].Engineering Mechanics，2008，25（10）：54－60.（In Chinese）

[11] KAREEM A.Numerical simulation of wind effects：aprobabilistic perspective[J].Journal of Wind Engineering and Industrial Aerodynamics，2008（96）：1472－1497.

[12] WANG T G，COTON F N.Prediction of the unsteady aerodynamic characteristics of horizontal－axis wind turbines including three－dimensional effects[J].Proceedings of the Institution of Mechanical Engineers Part A：Journal of Power and Energy，2000，214（A5）：385－400.