白澤升， 王孟鴻

(北京建筑大學 土木與交通工程學院， 北京 100044)

實際風荷載通常包含長周期部分和短周期部分，其中長周期部分通常以平均風速描述，短周期部分以脈動風來描述。對結(jié)構(gòu)物來說，由平均風速引起的作用力相當于靜力作用，不同部位通過體形系數(shù)來計算；脈動風速引起的作用力是動力性質(zhì)的，將引起結(jié)構(gòu)振動，通過風振系數(shù)來計算。荷載規(guī)范[1]規(guī)定：對于基本自振周期T1大于0.25 s的工程結(jié)構(gòu)，如房屋、屋蓋及各種高聳結(jié)構(gòu),以及對于高度大于30 m且高寬比大于1.5的高柔房屋，均應(yīng)考慮風壓脈動對結(jié)構(gòu)發(fā)生順風向風振的影響。規(guī)范對于復雜的空間網(wǎng)架結(jié)構(gòu)的風振系數(shù)沒有給出相應(yīng)的數(shù)據(jù)和計算方法。因此，復雜網(wǎng)架結(jié)構(gòu)的風振計算通常需要借助相關(guān)理論，通過計算分析確定[2]。

目前對復雜網(wǎng)架結(jié)構(gòu)進行風振響應(yīng)分析，主要采取對結(jié)構(gòu)的周圍風場進行模擬，以獲取結(jié)構(gòu)表面風壓數(shù)據(jù)，以進行后續(xù)的風振效應(yīng)分析；而風洞試驗一般作為剛性模型，可對結(jié)構(gòu)表面風壓分布進行測試，直接測試結(jié)構(gòu)的風振效應(yīng)較為困難。對結(jié)構(gòu)周圍風場的模擬方法主要有風洞試驗法和數(shù)值風洞法。風洞試驗法周期長、成本高，試驗數(shù)據(jù)需要經(jīng)過復雜的數(shù)據(jù)處理后導入分析軟件進行分析計算。數(shù)值風洞法成本低、周期短、效率高，便于參數(shù)化控制，但需要復雜的前處理建模工作，而網(wǎng)架結(jié)構(gòu)形式日益多樣化、異型化，無疑增加了建模工作量，同時計算流體動力學(CFD)軟件與結(jié)構(gòu)分析軟件脫節(jié)，網(wǎng)架結(jié)構(gòu)桿件與節(jié)點眾多，需要大量的荷載施加工作。

本接口軟件可將復雜網(wǎng)架模型自動導入前處理軟件ICEM，節(jié)省了大量的建模工作，而后處理功能基于FLUENT計算風工程功能，自動完成數(shù)據(jù)的讀取以及荷載施加與計算過程，減少了大量手動施加工作。

1 風振分析方法

由于風荷載的隨機特性，目前大跨空間網(wǎng)架結(jié)構(gòu)的風振分析方法主要有頻域法和時域法。頻域法是通過Fourier變換將風壓時程轉(zhuǎn)化為風壓功率譜，再將風壓功率譜通過動力系數(shù)傳遞轉(zhuǎn)化為動力反應(yīng)譜，并根據(jù)隨機振動理論通過對動力反應(yīng)譜進行積分求解得到結(jié)構(gòu)的動力響應(yīng)。因頻域法計算效率較高，是目前大多數(shù)學者常采用的方法[3-6]，其進行的是線性范圍內(nèi)的分析。而大跨網(wǎng)架結(jié)構(gòu)具有頻率密集性，需通過振型分解解耦處理，考慮振型模態(tài)對系統(tǒng)應(yīng)變能的貢獻[8]，選取主要貢獻模態(tài)，可能會損失一部分精度。時域法則是將結(jié)構(gòu)風壓荷載時程作用于大跨網(wǎng)架結(jié)構(gòu)進行風振響應(yīng)動力時程分析，獲得動力時程響應(yīng)的全過程。由于時域法能夠反映風荷載的時空特性，并能夠用于結(jié)構(gòu)的非線性分析。因此，本文采用了時域法。

進行時域內(nèi)的風振分析，需要確定作用在單元節(jié)點上的風荷載時程，通常有2種方法。其一采用風洞試驗的風壓時程，也是本文所采用的方法：從FLUENT模擬中提取風荷載時程，并將其等效為節(jié)點動力時程作用在相應(yīng)的網(wǎng)架節(jié)點上。其二是對脈動風進行模擬，使模擬出的風速盡可能接近和滿足自然風特性如平均值、與高度有關(guān)的自功率譜和互功率譜以及相位角關(guān)系等，常用的方法有諧波疊加法和線性濾波器法。通過模擬得到的風速時程在準定常假定的基礎(chǔ)上轉(zhuǎn)換成風荷載的時程。如何模擬合成輸入風速以獲得合理的風壓時程文件需要做進一步的研究工作。

2 STADS風振分析流程

2.1 湍流數(shù)值模擬方法

CFD數(shù)值模擬方法可分為直接數(shù)值模擬法(DNS)、雷諾平均法(RANS)和大渦模擬法(LES)。直接數(shù)值模擬法(DNS)是直接采用瞬時的N-S方程對各尺度的渦流進行求解,需要耗費較高的計算資源，目前無法應(yīng)用于真正意義上的工程計算。雷諾平均法(RANS)將湍流運動時均化處理，引入各湍流模型對流體控制方程進行封閉，廣泛應(yīng)用于穩(wěn)態(tài)風場求解中。大渦模擬法(LES)的原理是將湍流中的渦旋分成大渦和小渦，用瞬時N-S方程對大渦直接求解，采用亞格子模型考慮小渦的能量輸運與耗散，實質(zhì)是對湍流運動進行空間平均化。由于大渦模擬法(LES)對高雷諾數(shù)的復雜湍流運動具有較高的分辨率，可以捕捉到復雜的非定常流動現(xiàn)象及其各物理量，準確度較高，是目前求解非定常流較為理想的湍流數(shù)值模擬方法。故本文采用了大渦模擬法(LES)。

2.2 瞬態(tài)場風壓計算

基于大渦模擬法(LES)的瞬態(tài)場風壓計算程序：

1) 模型轉(zhuǎn)換，將STADS的模型轉(zhuǎn)換到ANSYS ICEM CFD[8-9]。

2) 瞬態(tài)風壓計算，通過FLUENT計算并輸出控制點的風壓時程[8,10]。

3) 計算節(jié)點的風壓力時程，通過節(jié)點的控制面積和風壓計算。

4) 求解振型模態(tài)[11]計算結(jié)構(gòu)自振周期，利用前兩階振型計算阻尼系數(shù)。

5) 利用Newmark法進行時程計算，獲得各節(jié)點位移時程。

6) 根據(jù)上一步的位移時程，計算各個節(jié)點的風振系數(shù)。

3 算例分析

本文以跨度93 m、長度94 m的大型扁筒殼為例展示STADS風振計算的全部過程。在此需要說明：由于監(jiān)控點太多，造成計算時間過長。因此，除進行模型分區(qū)處理以外，程序還增加了計算輸出監(jiān)控點時程的歸并處理。通過定義數(shù)個相鄰分區(qū)共用一個區(qū)域的風壓時程[12]，減少數(shù)據(jù)輸出，從而加快程序運行。

計算模型如圖1所示：跨度93 m，長度94 m，節(jié)點數(shù)1 472，桿件5 724，上弦塊數(shù)為712；經(jīng)分區(qū)，上弦劃分為142塊區(qū)域如圖1所示。

圖1 計算模型Fig.1 Computational model

3.1 模型轉(zhuǎn)換與瞬態(tài)風壓計算

通過選擇風振控制項自動選擇大渦模擬法(LES)進行計算，根據(jù)ICEM的腳本語言和函數(shù)控制格式生成模型文件，根據(jù)FLUENT二次開發(fā)方式，生成UDF和腳本控制文件。相應(yīng)參數(shù)控制，如圖2所示。

圖2 參數(shù)控制Fig.2 Parameter control

其中，大渦模擬亞格子模型采用壁面適應(yīng)局部渦黏模型(WALE)，對流場的整體流動以及分離、附著、轉(zhuǎn)捩等局部流動模擬較為準確。壓力與速度耦合采用SIMPLEC算法，對流項采用有限中心差分格式離散，時間項離散采用有限二階精度格式，壓力方程采用二階精度離散格式。

瞬態(tài)分析中，時間步長的設(shè)置對大渦模擬的結(jié)果影響很大，決定了模擬結(jié)果的時間分辨度，是大渦模擬收斂情況的一個重要參數(shù)，同時也要保證數(shù)值模擬能在一個合理時間內(nèi)完成。針對以上標準進行試算，計算步長控制為0.06 s，滿足模擬需要。

在ICEM CFD中讀取模型文件自動完成模型的轉(zhuǎn)換和網(wǎng)格劃分等前處理功能，程序通過UDF定義邊界條件，運行FLUENT腳本控制文件后，自動完成模型的選擇、監(jiān)控點的設(shè)置、求解參數(shù)設(shè)置等操作，并輸出控制點的風壓時程數(shù)據(jù)。控制點風壓時程如圖3所示。100、137號為山墻端部節(jié)點，30、34、70號分別是筒殼迎風面、頂部和背風面節(jié)點。

圖3 風壓- 時間曲線Fig.3 Wind pressure-time curve

3.2 節(jié)點風壓力時程的生成

時程分析前需要將FLUENT計算得到的控制點的風壓時程轉(zhuǎn)換成網(wǎng)架節(jié)點的風壓力時程，該過程涉及復雜每個網(wǎng)格承載面積的計算，通過STADS-2加載程序自動完成該功能，從而形成每一個節(jié)點的風壓時程?？缰?52、769、770號節(jié)點風壓時程如圖4所示。

圖4 節(jié)點風壓- 時間曲線Fig.4 Node wind pressure-time curve

3.3 自振周期和阻尼系數(shù)計算

利用子空間迭代求得前兩階振型周期：T1=0.641 s，T2=0.529 s；何艷麗等基于各基本參數(shù)對某筒殼風振響應(yīng)進行了研究，隨著阻尼比增大，各節(jié)點風振系數(shù)逐步遞減[13]，參考其阻尼比對風振系數(shù)的影響和鋼結(jié)構(gòu)阻尼比取值經(jīng)驗，阻尼比取值ξ=0.025；阻尼矩陣采用Rayleigh阻尼矩陣：[C]=α*[M]+β*[K]，其中，[C]、[M]、[K]分別為阻尼、質(zhì)量和剛度矩陣，α和β為比例系數(shù)。程序據(jù)依前兩階角頻率算得：α=0.268 2，β=0.002 3，如圖5所示。

圖5 計算系數(shù)求解設(shè)置Fig.5 Calculation coefficient solution setting

3.4 時程分析計算

時程分析計算采用高級分析軟件STADS-3進行。程序讀取節(jié)點風荷載時程，采用Newmark法進行時程分析，計算得到各節(jié)點位移時程，以用于下一步計算各個節(jié)點的風振系數(shù)。計算所得的跨中727和778號節(jié)點z向位移時程如圖6所示。

圖6 節(jié)點位移- 時間曲線Fig.6 Node displacement-time curve

3.5 風振系數(shù)的計算

根據(jù)荷載規(guī)范[1]定義，在工程應(yīng)用中將風荷載的動力效應(yīng)以風振系數(shù)β的形式等效為靜力荷載，即風振系數(shù)為：

(1)

采用峰值因子法，由平衡方程，荷載風振系數(shù)通過剛度矩陣和位移表示為：

(2)

對于節(jié)點i，荷載風振系數(shù)βi為：

(3)

式中：j為i節(jié)點的相鄰節(jié)點，kij為總剛矩陣第i行的非對角元。

采用峰值因子法，可將節(jié)點i處的位移風振系數(shù)表示為：

(4)

上述計算過程通過菜單引導執(zhí)行，將各節(jié)點風振系數(shù)記錄在相應(yīng)文件中。用戶可以選擇執(zhí)行荷載風振系數(shù)式(3)或者位移風振系數(shù)式(4)，其中在計算荷載風振系數(shù)時需要調(diào)用整體剛度矩陣。筒殼表面和山墻的位移風振系數(shù)如圖7所示。

3.6 風振系數(shù)的利用

上述無論荷載風振系數(shù)還是位移風振系數(shù)，都僅作為設(shè)計參考，需要設(shè)計人員做出合理判斷使用?？紤]到其加載過程手工施加的工作量較大，在STADS-2的加載過程中可以自動讀入各節(jié)點風振系數(shù)數(shù)據(jù)，可避免手工輸入的煩瑣與可能產(chǎn)生的數(shù)據(jù)錯誤，減小了相應(yīng)工作量。

4 結(jié) 論

1) 在STADS網(wǎng)架設(shè)計軟件基礎(chǔ)上，借助FLUENT軟件的風工程計算功能，編制了軟件模型轉(zhuǎn)換和風壓讀取接口，并在高級分析軟件STADS-3的動力時程分析基礎(chǔ)上，完成了復雜體型網(wǎng)架結(jié)構(gòu)的風振分析。本接口軟件的開發(fā)可為此類鋼結(jié)構(gòu)軟件的風振分析開發(fā)提供理論依據(jù)與思路參考。

2) 接口程序通過菜單引導控制相應(yīng)參數(shù)自動完成模型轉(zhuǎn)換、FLUENT風場數(shù)值模擬參數(shù)設(shè)置和大量控制點數(shù)據(jù)的輸入與計算結(jié)果的提取，以及荷載的自動讀取與施加，操作簡便，減少了大量手動輸入的煩瑣工作，節(jié)約了數(shù)值風洞環(huán)節(jié)的時間和經(jīng)濟成本。將FLUENT與STADS建立接口，有效彌補了計算流體動力學軟件(CFD)與結(jié)構(gòu)分析軟件的脫節(jié)，加快了程序的運行，提高了工作效率。

3) 接口程序以某實際工程的筒殼模型為算例，計算得到網(wǎng)架結(jié)構(gòu)各節(jié)點位移風振系數(shù)或荷載風振系數(shù)，在實際工程設(shè)計中可通過靜力加載方式計算得到結(jié)構(gòu)響應(yīng)，簡化了工程設(shè)計工作，為復雜體型網(wǎng)架結(jié)構(gòu)抗風設(shè)計提供參考。同時可適應(yīng)各類異型空間網(wǎng)架結(jié)構(gòu)，具有良好的實用前景。