考慮材料非線性的RC 雙曲冷卻塔風致破壞過程

張軍鋒，裴 昊，朱 冰，劉慶帥

(鄭州大學土木工程學院，鄭州 450001)

RC 雙曲冷卻塔作為一種高聳空間薄壁結構，風荷載是其設計控制荷載[1]。渡橋電廠冷卻塔風毀事故[2]和不斷攀升的結構高度使結構的抗風安全性始終備受關注，既有研究主要關注結構的風荷載特性、群塔干擾效應、風振效應以及彈性穩定問題[3]。實際上，渡橋電廠事故被一些學者認為是風致失穩破壞，從而引發了近四十年的穩定研究熱潮[4 ? 8]。然而，Mang 等[9]于1983 年采用有限元方法對Port Gibson 冷卻塔首次進行了考慮材料非線性的極限荷載分析，并認為其風致破壞應屬于材料破壞而非彈性失穩。之后，冷卻塔的靜風極限荷載分析受到關注。既有研究基本均對塔筒采用分層殼單元，沿厚度劃分為多個層次并用鋼板層模擬鋼筋。對于荷載主要考慮自重和風荷載[9 ? 13]，部分學者亦考慮了溫度效應[14 ? 15]。但由于問題本身的復雜性以及各學者所采用的軟件、單元、材料本構、計算方法以及荷載取值等的差異，即使對Port Gibson 同一冷卻塔的計算結果也表現出較大的離散性[16]：所得開裂荷載系數和極限荷載系數分別分布在0.75～1.58 和1.49～3.15，并且其中亦有文獻[11]認為是失穩破壞。

近年來，隨著國內冷卻塔建設規模的快速發展，國內學者亦注意到對冷卻塔應考慮材料非線性進行分析。但初期研究往往比較粗略，未能對混凝土的材料非線性進行準確模擬[17 ? 18]。文獻[19]進行了考慮初始缺陷的極限風荷載分析，但未就關鍵的材料非線性模型和模擬方法進行介紹。另外，對塔筒破壞的原因也有不同的結論：文獻[20]關注了動力極限風荷載，但僅對塔筒的受壓損傷進行分析而對更關鍵的受拉損傷卻沒有說明，并認為迎風區混凝土的破壞是受壓破壞；文獻[21]模擬了冷卻塔的倒塌過程，其結果顯示冷卻塔破壞源于混凝土受拉開裂和環向鋼筋斷裂引起的材料強度失效，這均與大多認為冷卻塔源于迎風區子午向受拉破壞的結論[9 ? 10, 12]不一致。鑒于塔筒RC 材料非線性模擬的困難，文獻[22 ? 23]將塔筒視為均質材料并賦予彈塑性材料特性，但這種處理無法準確模擬混凝土開裂后的特性以及鋼筋和混凝土之間的協同工作特點，且兩者對同一冷卻塔所取的彈塑性參數也差異顯著。

以一座大型冷卻塔為例，采用ABAQUS 軟件對其進行極限靜風荷載計算。采用分層殼單元模擬塔筒，采用彌撒開裂模型模擬混凝土的非線性受壓特性、受拉開裂和拉伸硬化效應，并計入幾何非線性。從荷載-位移曲線、變形模式、內力、應變和應力的發展對塔筒的破壞過程進行系統分析，并從中探究其破壞機理。

1 材料非線性模型

1.1 混凝土彌散開裂模型

鋼筋混凝土的承載特性表現出高度非線性，合理地反映這些非線性因素是準確獲得冷卻塔破壞過程的必要前提。本文采用彌散開裂模型模擬混凝土非線性本構，該模型可以準確地模擬混凝土構件受彎時的非線性受力特性[24]，也是既有冷卻塔極限風荷載分析的主要方式[9 ? 15]。

該模型通過4 個方面的材料行為模擬混凝土的非線性特性：非線性受壓特性，破壞面，剪切保留和拉伸硬化，詳見文獻[16]。

剪切保留模型如式(1)[25]和圖1[13]所示。

式中：G0為閉合混凝土的剪切模量；G'為開裂混凝土的剪切模量；β 為裂縫閉合后的剪切剛度系數；εt、εcr和εm分別為混凝土開裂應變、開裂后垂直于裂縫方向的應變和剪切剛度為0 時的應變。參考文獻[12 ? 13]，下文取β=0.25，εm=0.004。

拉伸硬化效應[10 ? 12]用以模擬混凝土開裂后的行為[24]。與文獻[10 ? 11]一致，下文同樣假定混凝土開裂后應力線性卸載(圖2)，在應變達到極限應變εu=TSεt時應力降至0，TS 即為拉伸硬化參數。既有研究[10 ? 12]表明TS 對結果影響較大，在板殼結構試驗結果驗證模型的基礎上[16]，下文取TS=15。文獻[16]還分析了TS 對結果的影響。

圖1 剪切保留模型Fig.1 Shear retention model

圖2 拉伸硬化模型Fig.2 Tension stiffening model

1.2 鋼筋非線性模型

盡管鋼筋被模擬為平面應力鋼板層，但其實際是單向受力，故采用單軸雙線性應力-應變曲線模擬其非線性本構[9,12,15]，如圖3 所示。其中：σy和εy分別為屈服應力和屈服應變；Es和Ep分別為彈性模量和塑性模量，下文取Ep=0.05Es。

圖3 鋼筋雙線性模型Fig.3 Bilinear model for reinforcement

2 結構及有限元模型

本文所選冷卻塔如圖4 所示，塔高177 m。下文采用塔筒相對高度hs/Hs表示所關注的高度位置hs，其中Hs為塔筒總高度(圖4)。塔筒為漸變厚度(圖5(a))，最小厚度0.271 m，塔筒下緣最大厚度1.40 m，上緣厚度漸增為0.40 m，塔筒通過48 對Φ1.30 m 人字柱與基礎連接。塔筒和頂端檐口均采用C40 混凝土和HRB400 級鋼筋，塔筒環向和子午向配筋率如圖5(b)所示，下支柱采用C45 混凝土。

圖4 冷卻塔特征尺度示意Fig.4 Characteristic dimensions of cooling tower

圖5 塔筒結構參數Fig.5 Structural parameters of tower shell

依規范[26]取基本風速V0=31 m/s，場地類別為B 類，對應風振系數β=1.9。溫度和風荷載聯合作用下的環向彎矩亦需重視[27]，考慮到冬溫和夏溫效應接近[28]，故僅考慮冬溫作用：依規范[26]塔外設計氣溫取?9.6 ℃，塔內設計氣溫取(環梁有擋水設施+單元系統)模式。

由于結構和荷載均對稱，為降低計算代價，僅建立半結構模型[9? 12,14]。塔筒采用殼單元S4R模擬，經不同網格密度劃分對比將塔筒劃分為環向×子午向=96×72 個單元，并沿厚度劃分9 個混凝土層加4 個單向受力鋼板層(圖6)。

圖6 分層殼示意圖Fig.6 Layered shell element

其中，最外層(也即圖6 中的第1 和第9 層)厚30 mm 用以模擬混凝土保護層，其它各層按1∶1∶1.5∶2∶1.5∶1∶1 確定。不同網格和層次劃分的對比分析詳見文獻[16]。鋼板層厚度根據配筋面積相等進行確定，其位置均在相應鋼筋的中心位置。下支柱和頂端檐口均采用梁單元B31模擬。

考慮兩種荷載工況并采用牛頓迭代法進行計算：自重+風荷載(G+λW)和自重+冬溫+風荷載(G+T+λW)，其中G、T 和W 分別為規范[26]所列荷載標準值。對于荷載效應，軸力F 以受拉為正，彎矩M 以外表面受拉為正；軸力和彎矩均以X 方向為環向，Y 方向為子午向，如FY和MX分別為子午向軸力和環向彎矩，具體可參考文獻[30]。徑向位移以內凹為正，外凸為負。環向角度θ 以迎風點為0°，逆時針轉動為正。

需要說明，在對本算例進行分析前，還首先根據文獻[10, 12]對Port Gibson 冷卻塔進行了重復計算[16]。盡管所用軟件和單元不同，材料非線性模型和計算方法等也存在差異，但所得結果基本吻合(圖7)，從而驗證了本文所用軟件、單元以及材料非線性模型和計算方法的合理性。

圖7 Port Gibson 塔的荷載位移曲線Fig.7 Load displacement curves of Port Gibson tower

3 線性結果分析

作為非線性分析的基礎，首先應了解冷卻塔在各荷載作用下的內力特征。自重作用下(圖8)塔筒以雙向受壓為主，且子午向壓力遠大于環向壓力，雙向彎矩可以忽略。自重產生的軸壓力將有助于抵抗風荷載產生的拉力。冬溫荷載作用下(圖9)塔筒以雙向正彎矩為主，雙向軸力可以忽略。

圖8 自重作用下塔筒內力分布Fig.8 Internal forces of shell under gravity

圖9 冬溫荷載作用下塔筒內力分布Fig.9 Internal forces of tower shell under gravity

風荷載作用下塔筒內力分布較為復雜(圖10)，除MY外，各荷載效應均較顯著，并以子午向軸力最為顯著。對各內力特征的詳細介紹可見文獻[27]。

圖10 風荷載作用下塔筒內力分布Fig.10 Internal forces of tower shell under wind load

4 非線性結果分析

4.1 變形模式

在兩種工況下，施加風荷載前，結構變形可以忽略。施加風荷載后，各截面的變形規律基本一致，并且考慮到對稱性僅給出0°≤θ≤180°范圍的結果(圖11)：即迎風區和側風區分別產生顯著的內凹和外凸位移，且前者比后者增加更為顯著，背風區的變形始終較小。在進入非線性階段后，迎風區內凹位移最大值均始終在θ=0°位置，但隨風荷載的增加而有上下移動；而側風區外凸位移最大點則同時有環向和子午向移動：圖11 所示從下到上的4 個斷面的外凸位移最大點逐漸由θ=40°區域后移至θ=75°區域，而在子午向則逐漸由hs/Hs=0.65 區域上移至塔筒頂部。

圖11 塔筒橫向截面變形Fig.11 Deformed shapes of tower shell: sliced views

4.2 荷載-位移曲線

為直觀體現結構的非線性效應，選?、顸c(θ=0°、hs/Hs=0.782，迎風點喉部)和Ⅱ點(θ=60°、hs/Hs=0.70)這兩個內凹和外凸位移始終較為顯著的位置，以及Ⅲ點(θ=120°、hs/Hs=0.70)作為對照位置，給出荷載-位移曲線如圖12 所示?？梢钥闯?，兩種荷載工況下，冷卻塔破壞過程中的非線性主要體現在迎風區和側風區，其他區域的非線性極不明顯。

圖12 兩種荷載工況所得荷載位-移曲線Fig.12 Load-displacement curves of two load cases

對于工況G+λW，以Ⅰ點為例，其荷載位移曲線可劃分為3 個階段。第1 階段(OA 段，λ<1.384)為線性響應階段：位移隨風荷載線性增加，整個塔筒未出現裂縫。第2 階段(AB 段，1.384≤λ<1.876)為裂縫擴展階段：由于風荷載在塔筒迎風區產生的子午向軸拉力逐漸達到并超過自重作用下的子午向軸壓力，故λ=1.384 時，迎風子午線上hs/Hs=0.37 位置混凝土應力首先達到峰值并進入卸載階段，也即該處產生裂縫，且該裂縫為貫穿塔筒壁厚的環向裂縫；隨著風荷載的增加，環向裂縫在塔筒迎風區不斷形成并擴展，塔筒局部剛度逐漸下降，結構響應進入明顯的非線性階段，迎風區和側風區位移迅速增大，內力也出現重分布，詳見下節。第3 階段(BC 段，1.876≤λ≤2.007)為硬化階段：λ=1.876 時，迎風子午線上hs/Hs=0.70位置內側子午向鋼筋首先受拉屈服，結構進入硬化階段，風荷載又有進一步增加，鋼筋屈服區域沿環向逐漸擴展至θ=16°，沿子午向則擴展至hs/Hs=0.52～0.80。并且λ=1.529 時，迎風子午線hs/Hs=0.62高度內側和θ=75°、hs/Hs=0.95 高度外側出現子午向裂縫。最終，由于迎風區混凝土的持續開裂和鋼筋屈服，λ=2.007 時計算終止(表1)。這一荷載位移曲線形式與文獻[9 ? 15]所得結果較為類似。

表1 兩種荷載工況結果Table 1 Results of two load cases

對于工況G+T+λW，由于冬溫荷載產生雙向正彎矩，使塔筒外表面混凝土在風荷載作用前即處于雙向受拉狀態，而風荷載作用后側風區的環向正彎矩效應進一步增加，迎風區子午向也將處于拉彎聯合受力狀態，所以塔筒提前開裂：λ=1.000 時，塔筒側風區θ=73°、hs/Hs=0.93 位置外表面即首先出現子午向裂縫；λ=1.080 時，塔筒迎風子午線上hs/Hs=0.08 位置外表面開始產生環向裂縫。因此，該工況的荷載-位移曲線略提前進入非線性階段，且整體剛度偏小。

對比表1 可知，冬溫荷載顯著減小了冷卻塔的開裂荷載和極限位移，但對極限荷載影響不大。這是由于盡管冬溫荷載產生的雙向正彎矩導致塔筒提前開裂，但此彎矩效應較小且固定不變，且隨著風荷載的持續增加，其與風荷載產生的荷載效應相比也越來越小，故冷卻塔的最終破壞依然由風荷載控制。也正因此，兩種荷載工況所得荷載-位移曲線幾乎一致；二者變形和裂縫開展情況也相似(見4.2 節和4.3 節)，只是在同一λ 下，計入冬溫效應的變形和裂縫開展更為嚴重。因此，對于本例，溫度作用只是增加了環向和子午向彎曲損傷，這在文獻[14 ? 15]中亦有體現。

4.3 內力分布

荷載-位移曲線表現出的明顯非線性特征直接源于混凝土開裂導致的結構剛度損傷，這也使結構剛度分布發生變化，進而導致內力重分布并加劇結構開裂。由于兩種荷載工況所得荷載-位移曲線變形模式相差不大，為觀察塔筒破壞過程中的內力重分布情況，以工況G+λW 為例，給出hs/Hs=0.37 和0.70 這兩個斷面的內力分布如圖13 所示。

由圖13 可知，在λ≤1.384 時，整個塔筒所有內力均呈線性變化，并且由于計入了自重效應，FY在迎風區逐漸由軸壓力表現為顯著的軸拉力，在側風區始終為顯著的軸壓力，在背風區亦為軸壓力但幅值較??；自重效應的存在也使FX在整個環向均為軸壓力；MX在整個環向幅值較小，在迎風區和側風區分別為負彎矩和正彎矩。

圖13 塔筒內力分布Fig.13 Distribution of internal forces of tower shell

λ=1.384 時，由于hs/Hs=0.37 高度首先出現環向裂縫(見下文4.5 節)，θ=0～15°的各內力均開始表現出非線性并逐漸擴展到θ=45°位置，即出現內力重分布。并且環向裂縫的出現不僅影響子午向受力，也影響環向受力：與此同時，hs/Hs=0.70 高度θ=0°的FX和MX也開始出現內力重分布，而該高度FY依然保持線性變化，直至λ=1.48 時開始非線性變化。而此時此處尚未開裂，這同樣是內力重分布的結果。此后，隨著塔筒開裂的繼續發展，塔筒內力重分布顯著區域也不斷擴展至θ=90°位置并在兩個高度表現出相似的規律。

對于子午向的內力重分布，θ=0°位置開裂后，此處的FY略有下降，但其兩側的FY快速增加，用以補償θ=0°子午線上FY的下降。迎風區大范圍受拉開裂后，整個迎風區的FY不再線性增加，還造成側風區的FY亦不再隨荷載線性增加。與環向裂縫主要因子午向軸拉力引起不同，子午向裂縫主要由環向彎矩引起，內力重分布也主要表現為迎風區和側風區MX和FX幅值和環向波動的急劇增加。

4.4 應變分布

圖14 給出了兩種荷載工況所得塔筒最外層和最內層混凝土最大主應變εpr分布。由于冬溫和風荷載的彎矩效應有限，兩種工況下，塔筒均主要依靠薄膜作用抵抗風荷載，所以塔筒內外層混凝土εpr隨風荷載的變化規律基本相同，尤其是工況G+λW。對于工況G+λW，λ=1.50 時，迎風區中下部已有較大面積混凝土產生了環向裂縫，故結構開始表現出明顯的非線性特征。λ=1.80 時，迎風區已產生大面積開裂，并存在4 個明顯的貫穿壁厚的環向裂縫帶(亦可見圖15(c))；由于迎風區和側風區分別存在較大的環向負彎矩和正彎矩，故上述區域內側和外側混凝土均產生了子午向裂縫，但其寬度相對較小，直到結構破壞也未貫穿壁厚。λ=2.00 時，迎風區hs/Hs=0.08～0.80 區域全部開裂，且開裂區域與λ=1.90 時近乎一致，僅裂縫寬度急速增大。最終，冷卻塔因大面積的混凝土開裂和鋼筋屈服而破壞。計入冬溫荷載后，冬溫荷載產生的雙向正彎矩導致塔筒提前開裂，且同一λ 下裂縫開展情況更為嚴重，主要表現在側風區外側。

4.5 應力分布

由于迎風子午線上的子午向受拉開裂和hs/Hs=0.70 斷面的環向拉彎開裂最為顯著，為進一步了解塔筒破壞過程中的應力分布情況和破壞過程，給出G+λW 工況下迎風子午線上最外層混凝土和外側鋼筋的子午向應力分布和hs/Hs=0.70 斷面表層混凝土和鋼筋的環向應力分布，同時給出上述兩個位置處的應力荷載曲線，如圖15 所示。

由圖15(a)可以看出，λ=0.50 時，除塔筒頂部和底部邊界外，塔筒整個高度范圍均處于雙向受壓狀態。λ=1.0 時，風荷載作用下的子午向軸拉力已經超過自重的壓力效應，迎風區大部分區域已進入子午向受拉狀態。λ=1.50 時，hs/Hs=0.10～0.50范圍內混凝土子午向應力已達到抗拉強度，且hs/Hs=0.37 處混凝土已進入卸載階段。λ=1.80 時，迎風區已有4 個位置混凝土應力完全降至0，即迎風區產生了4 個完全開裂的環向裂縫帶，裂縫區荷載完全由鋼筋承擔導致鋼筋應力激增(圖15(b))。λ=2.00 時，迎風區子午向開裂區域幾乎不再隨λ增加，hs/Hs=0.08～0.80 范圍混凝土完全退出工作，hs/Hs=0.52～0.80 區域鋼筋屈服(圖15(b)、圖15(c))。

圖14 塔筒混凝土最大主應變(僅給出εt=7.354×10?5～5.2×10?3 范圍內的拉應變)Fig.14 Maximum principal strain of concrete in tower shell

圖15 塔筒外側混凝土和鋼筋應力分布Fig.15 Distribution of stress of exterior concrete and steel in tower shell

另外，從圖15(c)更可直觀看出，混凝土開裂后其應力急劇減小至0，與此同時裂縫區鋼筋的應力則急劇增加。4 個主裂縫帶產生的位置順序依次為hs/Hs=0.37、0.52、0.17 和0.70，對應的λ 依次為1.384、1.451、1.525 和1.631。盡管喉部位置混凝土最后開裂，但由于該處壁厚和配筋均較小，故該處鋼筋首先屈服。這也說明，子午向受拉開裂和鋼筋屈服是冷卻塔在風荷載作用下破壞的直接原因。

對于環向應力，由圖15(d)可以看出，hs/Hs=0.70高度外側混凝土環向應力的分布與圖13 所示MX的分布規律一致，這也表明冷卻塔的環向應力狀態受彎矩控制。λ≤1.50 時，除迎風子午線局部區域，環向應力基本處于線性狀態；進入非線性狀態后，受環向正彎矩影響，θ=25°～65°區域外側混凝土應力達到抗拉強度，但此后其應力并未快速卸載，而是維持在一個較高的水平，這是由于該區域的裂縫寬度在塔筒破壞過程中始終較小(圖14)，這也體現在該區域鋼筋應力在混凝土開裂后增長不大(圖15(e))。而由圖15(f)可知，λ=1.529和1.696 時，hs/Hs=0.70 高度θ=0°內側和θ=50°外側混凝土環向應力先后達到抗拉強度并進入卸載階段，這兩處混凝土產生子午向裂縫，且鋼筋的拉應力迅速增加。盡管這兩處環向鋼筋均未屈服，但λ=1.98 時迎風子午線上hs/Hs=0.57 處內側環向鋼筋屈服，并迅速擴展至hs/Hs=0.53～0.61 范圍。這也表明，盡管結構破壞源于子午向受拉破壞，但環向受力亦需關注。

5 結論

以一座大型RC 冷卻塔為例，通過數值模擬研究了冷卻塔在兩種荷載工況下的非線性靜風響應，從荷載位移曲線、變形模式、內力、應變和應力分布對其破壞過程進行了系統闡述，并從中探究了其破壞機理。主要結論如下：

(1)通過對Port Gibson 冷卻塔的計算和與既有結果對比，表明ABAQUS 的分層殼單元配合彌散開裂模型可以有效進行冷卻塔極限靜風荷載分析。

(2)在自重+風荷載作用下，在λ=1.384 時，迎風點hs/Hs=0.37 位置在子午向軸拉力作用下首先開裂并產生貫穿厚度的環向裂縫，之后裂縫隨風荷載增加而不斷形成并沿環向擴展至整個迎風區，沿子午向幾乎擴展至整個高度，開裂位置的鋼筋應力迅速增加；λ=1.529 時，迎風子午線hs/Hs=0.62高度內側和θ=75°、hs/Hs=0.95 高度外側混凝土因環向彎矩而開裂；λ=1.876 時，hs/Hs=0.70 處子午向鋼筋首先受拉屈服，之后鋼筋屈服區域沿環向逐漸擴展至θ=16°，沿子午向則擴展至hs/Hs=0.52～0.80；λ=1.98 時迎風子午線上hs/Hs=0.57 處內側環向鋼筋屈服，并迅速擴展至hs/Hs=0.53～0.61 范圍。λ=2.007 時，因混凝土持續開裂和鋼筋屈服，冷卻塔達到抗拉極限狀態而破壞。

(3)在自重+風荷載作用下，塔筒開裂之前各響應均表現為線性特征，塔筒開裂之后，塔筒局部剛度下降，荷載-位移曲線也由線性進入非線性并且迎風區和側風區位移迅速增加，同時迎風區和側風區內力也表現出明顯的重分布特征。結構破壞時，迎風區子午向受拉導致大面積的受拉貫通環向裂縫和鋼筋屈服，而環向僅有塔筒中上部受彎開裂且裂縫寬度有限，鋼筋屈服范圍也較小，因此可認為迎風子午向受拉導致的塔筒開裂和鋼筋屈服是結構破壞的首要原因，但環向受力亦需關注。

(4)在自重+冬溫荷載+風荷載作用下，由于冬溫荷載產生的雙向正彎矩效應增加了塔筒的彎曲損傷，導致塔筒側風區在λ=1.0 時提前產生子午向裂縫。但由于此彎矩效應固定不變且遠小于風荷載作用下的子午向受拉效應，故結構的破壞依然由風荷載控制，極限荷載系數為λ=1.842。

(5)兩種荷載工況下結構均在λ=1.50 時明顯進入非線性段，且之后塔筒變形急速增加，此荷載系數這也與結構設計時的風荷載分項系數1.4 相當，說明在設計風荷載下結構的冗余安全度較為有限。