基于流固耦合的骨架膜結構優化

李延民，黃 瑞，趙樹森，王景玉

(鄭州大學 機械與動力工程學院，鄭州 450001)

骨架膜結構是一種新型圍護結構，預緊后的膜材與金屬骨架整體構成一個高次超靜定空間結構，依靠膜材與骨架整體形成的結構剛度來承受外部荷載，在力學上屬于空間梁膜混合結構，具有典型的非線性特征.對比傳統圍護結構，具有結構簡單、自重低、造型多樣、抗震性能優良的特點，其制造規格統一，多為工廠化生產后現場拼裝，便于快速搭建一個大跨度無支撐結構，能提供更大的有效空間[1].廣泛應用于工業、建筑、國防等領域中作為設備、儀器、建筑的圍護結構[2].

骨架膜結構在服役生命周期中暴露于自然環境下，受多種環境荷載的影響，其中最為嚴重且最不規律的是風荷載[3].因此，有必要對此類圍護結構進行風場分析，目前，該類結構設計時多按照相應規范取單一風向角(0°風向角)下的風壓分布進行設計驗算，不能反映出實際的最惡劣工況.

劉奕等[4]采用風洞測壓試驗對不同地貌下不同深寬比的矩形高層建筑進行了風荷載及其周身繞流特性的研究，得到了矩形高層建筑表面的一些風壓分布規律；李志國等[5]基于風洞試驗對某機場航站樓的大跨度曲面屋蓋的脈動風荷載分布進行了研究，得出屋蓋迎風前緣及部分曲面弧度變化大的位置處風荷載在不同風向角下呈現非高斯特性.但由于風洞試驗費用高、周期長，較少對膜結構進行全尺寸的風洞試驗，而縮尺模型由于雷諾系數的不同，其風場往往與實際風場有差別.Amaya-Gallardo、Amini、方智遠等[6-8]采用數值模擬的方法，與實尺模型場地實測數據、風洞試驗數據對比，證明數值模擬方法得到的建筑物表面平均風壓能較好地與場地實測數據、風洞試驗數據吻合.

現階段，國內外雖存在一些學者對低矮類建筑結構或膜結構的風場特性進行研究，但多將建筑結構本身視為剛性體，忽略了其本身在風荷載作用下的結構變形對其表面風壓分布的影響.流固耦合的相關概念雖然提出得很早，但是其對計算機的性能要求很高，限于早期計算機的性能，基于雙向流固耦合的風場特性研究鮮有報道.目前，隨著科學技術的發展，采用弱流固耦合計算的門檻已經降低.本文基于ANSYS，考慮了膜結構在風荷載作用下的結構變形對建筑表面風壓分布的影響，采用雙向流固弱耦合，對低矮類膜建筑的風場特性與骨架膜結構風致結構受力進行了分析研究.同時，骨架膜結構的優化一般通過優化桿件的截面尺寸來實現輕量化設計，但未考慮其設計荷載與實際荷載之間的差異.本文基于流固耦合所得到的多個風向角下建筑表面的風壓分布數據，對該建筑的立面骨架膜結構進行了輕量優化，相比減重14.73%，有效節約了原材料.

1 雙向流固耦合分析流程

針對此類對風荷載較為敏感的結構，采用雙向流固耦合分析與實際工況更為接近，與單向耦合分析不同的是，在需要耦合的流場計算表面區域創建動網格.將CFD求解器解出的耦合面風壓分布與FEM求解器解出的結構位移量，通過System Coupling模塊在每個時間步間互相交換作為各自輸入變量再進行下一個時間步的求解，依次循環，實現雙向流固耦合，這樣大大減少了將流場N-S方程與動力學方程聯合求解的巨大工作量且保證了求解的準確性，其求解流程如圖1所示.

圖1 雙向流固耦合分析的一般流程

2 風場分析

2.1 分析模型建立

該骨架膜圍護結構安裝示意圖如圖2所示(單位：mm).該結構總長80 m、寬30 m、高12 m，屋蓋坡度為15°；屋面設有4面整體框架均布平面骨架膜結構，長35 m、寬10 m，立墻上設有16面框架均布骨架膜結構，長10 m、寬6 m.為了便于模塊化制造與包裝運輸，將骨架膜結構設計為均勻大小，其中屋面骨架膜模塊尺寸為3 m×2 m、立面骨架膜模塊尺寸為2 m×1.5 m.

圖2 骨架膜結構安裝示意圖

相關調查研究表明，低矮類建筑在強風下，屋面與一些圍護結構處會產生局部低氣壓，在內外壓強差的作用下導致建筑物的圍護結構受損[9-10].因此有必要對不同迎風角下建筑表面的風壓分布進行研究，確定最危險的迎風角度.

2.2 計算域與網格模型

由于實際湍流場的參數分布規律難以獲得，本文使用充分發展的流場去模擬真實入口邊界條件，為了保證流場的充分發展，需要確定一個足夠大的流場區域，在風向角為0°(與YOZ平面夾角)的情況下，流場計算區域大小為600 m×448 m×240 m，流場入口距離建筑180 m，出口距離建筑后邊緣390 m[11].圖3為流場網格模型.采用分區組合的方式離散流場，遠離建筑的外區域采用六面體結構化網格.考慮到湍流區域空氣流向變化劇烈，在內區域采用非結構化的四面體網格，并在建筑周圍進行局部加密.整個計算流場網格模型共1 300 257個節點、2 891 472個單元.

圖3 流場網格模型

2.3 流場邊界條件

在Fluent中將湍流模型設定為k-εRNG模型，將入口的均勻流與實際湍流場建立聯系.根據統計，上海市50年重現期基本風速為31 m/s，因此入口風速設置為31 m/s，出風口壓力為標準大氣壓[12].近壁面處理方式為Standard Wall Functions；材料假定為不可壓縮粘性空氣，質量密度取1.225 kg/m3，粘度為1.789 4×10-5kg/m·s.

對于一個充分發展的流場模型，其湍流特性可以采用湍流強度I與湍流尺度l描述，湍流強度主要靠雷諾數Re估算，而湍流尺度與流場中對流動影響最大的因素有關，此處為建筑來流方向特征尺度，表達式為

(1)

I=0.16(Re)-1/8

(2)

l=0.07L

(3)

式中：ρ為流體密度；μ為流體動力粘度；u為流速；L為建筑來流方向特征尺度.將式(1)、(2)聯立得到湍流強度I=1.69%，由式(3)得到湍流尺度l=2.1 m.

3 風場數值模擬結果與驗證

3.1 風場數值模擬結果

通過仿真得到的0°風向角下的風場壓力分布如圖4～6所示.圖4為距建筑側墻5 m處剖面的風速矢量圖，可見由于氣流流向存在截面大小的劇烈變化，在迎風面方向屋檐與屋脊位置出現了明顯的氣流分離現象，導致局部產生了高負壓，而背風面由于建筑阻擋形成負壓，引起周圍空氣的回流，形成渦流致使背風面呈現局部正風壓.由圖5可見，計算區域囊括了整個受建筑物阻礙所形成的湍流區域，說明流場的大小滿足計算需求.由圖6可見，整個建筑迎風面受到來流撞擊產生正壓作用，其余區域均為風負壓作用.

圖4 風速矢量圖

圖5 中軸線截面風壓分布

圖6 建筑表面風壓分布圖

為了方便對比，將仿真得到的風壓分布數據以不受建筑物影響的遠高處風壓作為參考值進行無量綱化，得到平均風壓系數表達式為

(4)

式中：Q為建筑表面任意位置的風壓；QH為遠高處的靜壓；VH為平均風速.

圖7為建筑物屋面風壓系數分布等高線圖.

圖7 風壓系數分布等高線圖

3.2 結果驗證

為了驗證仿真得到的風壓數據是否有效，選取同類型低矮建筑模型TTU(美國德州理工大學)實尺模型場地實測數據[13](見圖8，單位：mm)、同濟大學風洞試驗數據[14]和建筑結構載荷規范(GB 50009-2012)中給出的低矮雙坡屋面建筑(坡度為15°)的風壓分布體系系數參考值，與建筑物模型中軸線上分布的測點處測得的風壓系數相比較，得到的對比曲線如圖9所示，測點位置排序按圖中A→B→C→D→E的方向散布.

圖8 TTU模型及測點布置圖

圖9 建筑物中軸線分布測點風壓系數對比

由對比結果可知，數值模擬結果與實測數據、風洞試驗數據具有相似的趨勢，與國標參考值吻合度較高，由于本文的建筑模型與TTU、同濟大學所采用的模型在尺寸和屋面坡度不同，故在迎風面、屋檐、屋脊處的風壓系數分布有一定的不同，總體而言，所得風壓數據合理準確，能夠滿足工程應用的需要.

類似的，改變風向角為45°、90°并設置對應的計算域大小和湍流參數，其余邊界條件保持不變，得到了45°與90°風向角下該建筑的表面風壓分布，如圖10所示.

圖10 45°、90°風向角下建筑物表面風壓分布圖

4 骨架膜結構風致結構受力分析

4.1 骨架膜結構幾何模型

本文以立面骨架膜結構為例，介紹該結構的流固耦合仿真與優化過程.骨架膜結構整體由膜材、上弦桿、腹桿、下弦桿、節點支撐桿與連接節點構成，整體受力分析時將該膜材簡化為shell181膜單元，上弦桿件簡化為beam188梁單元，腹桿、下弦桿與節點支撐桿簡化為只受拉壓的二力桿link180單元，為方便敘述，文中使用“二力桿”代指這些桿件[15].簡化后的模型如圖11所示.

圖11 骨架膜結構簡化模型

上弦桿是型號為8/5的角鋼，厚5 mm，腹桿與下弦桿采用截面為φ50×10 mm的鋼管，節點支撐桿是截面為φ60×10 mm的鋼管，金屬骨架的材質均為Q235，膜材厚度為0.8 mm，材質是型號為G7的PTFE膜材.

4.2 邊界條件

在Mechanical模塊中將膜面單元使用APDL命令施加6 kN/m的預張力，打開大撓曲開關，施加重力荷載，將膜面設置為流固耦合面，將膜結構的固定約束按照實際情況在上弦邊桿施加固定約束，在下弦支撐節點處施加不動鉸支座.離散后的屋蓋骨架膜結構節點總數為45 270，單元總數為41 674；立面骨架膜結構節點總數為28 661，單元總數為27 005.

4.3 仿真結果

在0°、45°、90°風向角下對立面骨架膜結構進行流固耦合仿真.圖12為骨架膜結構在不同風向角下的變形圖.表1為立面骨架膜仿真結果.

表1 立面骨架膜仿真結果

圖12 立面骨架膜結構在不同風向角下的變形圖

4.4 結果分析

由仿真結果可知，立面骨架膜結構在0°與45°風向角下承受的風荷載為正風壓，而在90°風向角下承受的為負風壓；在90°風向角下膜面變形的應力、骨架最大撓曲、受拉的二力桿最大軸力與桿件所受最大拉應力均比0°與45°風向角工況下大，但受壓的二力桿軸力與桿件所受最大壓應力卻產生于0°風向角下.其中，骨架最大撓曲位置為正中間的下弦桿中心，最大拉桿與最大壓桿均為兩側中間位置鉸支座處的四根斜腹桿，桿件所受的最大組合拉壓應力分別處于與該斜腹桿相連的上弦桿的上下翼緣處.膜面的最大應力差別不大，均遠低于膜材的極限抗拉強度標準值160 MPa.

在不同風向角下，風場的小范圍高負壓區域可能導致骨架膜圍護結構局部產生較大的變形撓曲、軸力和應力.從骨架膜結構在不同風向下結構受力情況來看，最大受力桿件為與支座鉸接的斜腹桿，其次是下弦桿，而骨架中間部位的腹桿總體受力較小，因此，可以考慮對這三種桿件使用不同截面尺寸的鋼管制造以輕量化整個結構.

5 骨架膜結構優化

在結構設計中，常采用響應面法對結構的相關設計參數進行優化，通過對大量樣本點進行計算，擬合出設計變量與目標變量之間的響應函數，并在給定的約束條件下挑選出最符合目標函數的候選點，將候選點圓整成制造值得到最終優化參數.一般優化流程如圖13所示.

圖13 響應面法的一般流程

根據建筑結構載荷規范對骨架膜結構取不同承載力極限狀態組合進行對比，其中風荷載采用仿真得到的不同風向角下骨架膜結構表面壓力分布數據，其余荷載按規范取值.根據對比可知，對結構整體承載最不利的荷載工況為：1.2D(恒載)-1.4W(90°風向角下負壓風載)+1.0F(膜面預應力)-0.7T(降溫)的組合.

以該極限工況組合作為施加荷載，使用ANSYS Workbench中的優化設計模塊對該骨架膜結構的金屬骨架進行截面尺寸優化，輸入參數為桿件的截面尺寸，其值取常用的角鋼、圓鋼管截面尺寸，輸出參數為結構的總質量，約束條件為金屬骨架最大撓曲、二力桿的長細比、壓桿穩定性與桿件的強度，優化數學模型為

(5)

式中：P1、P2、P3為上弦桿角鋼截面的長、寬、厚；P5、P7、P9、P11為下弦桿、節點支撐桿、支座處斜腹桿、中間腹桿的外徑；P6、P8、P10、P12為對應桿件的內徑；P13為骨架整體最大撓曲；P14為二力桿中受壓桿最大軸力；P15、P16為桿件所受的最大壓應力與拉應力；優化目標是使得P21最小.

使用Optimal Space-Filling Design法生成樣本點，刪除桿件內徑大于外徑的樣本點后共有151個樣本點.計算每個樣本點并使用響應面法對輸入、輸出參數進行擬合，得到了3個候選點，取3個點設計參數的均值并圓整后作為最優解，對所得的4個設計點進行計算驗證，如圖14所示.可見以最優解的驗證點為基準，其余7個點的誤差值均在10%以內，證明擬合度較高.

圖14 擬合誤差對比圖

圖15為各設計變量對輸出變量的影響.從圖15中可以看出：1)上弦桿、下弦桿和中間腹桿的截面尺寸對金屬骨架剛度影響最大，且為正相關；2)上弦桿和支座處斜腹桿的截面尺寸對二力桿中受壓桿所受軸力影響最大且呈現負相關趨勢，而下弦桿以及節點支撐桿的截面尺寸與二力桿所受軸力呈現正相關；3)加長上弦桿角鋼的截面長度能有效降低上弦桿所受壓應力，加寬截面寬度則能有效降低其所受拉應力；4)對結構整體質量影響最大的是弦桿角鋼的截面寬度和中間腹桿的內徑.

圖15 設計變量對輸出變量的影響

立面骨架膜結構的各桿件尺寸與總質量如表2所示.經優化后結構總質量為1 277.3 kg，降低了14.73%.

表2 骨架膜結構優化前后參數對比

6 結 論

本文通過分析得出以下結論：

1)基于ASNSY平臺采用雙向流固耦合方式對某項目中的骨架膜圍護結構進行不同風向角下的定常風場仿真，得到了結構在不同風向角下的表面風壓分布.

2)將所得風壓數據與TTU實尺模型場地實測數據、同濟大學風洞試驗數據和建筑結構載荷規范進行比照，驗證了數據的合理性與準確性.

3)采用雙向流固耦合方式將所得風壓耦合到骨架膜結構上，分析不同風向角下該建筑骨架膜結構各組成部件所受荷載情況，以此為據，將主要受力桿件與受力較小桿件的截面尺寸分開進行設計優化.

4)基于ANSYS的優化設計模塊，以不同風向角下的骨架膜結構表面風壓分布為基礎，選取對結構最不利的極限荷載工況組合作為施加荷載，將金屬骨架整體最大撓曲、二力桿的長細比、壓桿穩定性與桿件的強度作為約束條件，以最小化結構總質量為優化目標，對各類桿件的截面尺寸進行優化.

5)本文基于流固耦合技術與響應面優化法，得到了該骨架膜結構中各桿件的截面尺寸對結構受力的影響關系，并將該結構質量降低了14.73%.