深圳市陰香類行道樹風致失效模式及失效風速預測研究

李利孝，蔡 鳴，宋 健，李偉文，2，王衛侖，2

（1.深圳大學土木與交通工程學院，廣東 深圳 518060；2.深圳大學結構工程研究所，廣東 深圳 518060）

引言

據臺風年鑒統計，每年由臺風給我國東南沿海地區造成的經濟損失數以百億計，其中城市行道樹的破壞及其附帶損失（砸傷行人、車輛，飛射物對建構筑物的破壞等）是造成巨大經濟損失和次生災害的重要因素［1］。隨著人民對居住環境需求的日益提升，城市綠化生態建設持續增加，如深圳市的森林覆蓋率達到40%以上，綠化覆蓋率則超過50%。與此同時，臺風給深圳市行道樹造成的損毀也日益漸增，如2018 年第22 號臺風山竹粗略估計造成11 680 棵樹木倒損，經濟損失達數10 億。因此，對城市行道樹抗風性能開展系統研究，建立其風致破壞模式、確定其破壞風速及臺風易損性對城市綠化以及防臺減災具有重要指導意義。

國際上對于樹木抗風性能開展了一個多世紀的探索，早在1881年，Greenhill［2］將樹木簡化成勻質的錐形體，研究了樹木的風致振動特性。考慮到樹冠的自重，Spatz等［3］對Greenhill的模型進行了改進，在頂端附加了集中質量。Baker［4］則建立了樹木風致倒伏的簡化理論模型。Sellier 等［5］建立了樹木的有限元模型，推動了樹木風振響應研究的發展。Ciftci 等［6］使用有限元法討論了楓樹樹冠的空間形態對其動力放大因子的影響。Dupont等［7］建立了冠層的隨機風速模型，預測了不同風速下樹木的破壞概率。

在國內，對樹木風致破壞的研究主要集中在臺風過后的樹木受損情況調查，提出樹木倒伏的影響因素及減災措施。對于樹木抗風性能方面的研究，宮曉芳［8］利用編程結合分形理論建立樹木力學結構模型。艾曉秋等［9］采用頻域分析方法對行道樹的失效風速進行了研究，建立了考慮樹木形態和分枝影響的行道樹風災破壞評估有限元模型。李正農［10］采用多松弛格子Boltzmann 模型對樹木流場進行分析，準確預測了樹木最大風速的位置；隨后考慮了樹木間干擾作用，通過風洞試驗初步探明樹木間的風荷載與流場特性［11］。

本研究將旨在建立考慮樹木形態和力學性能影響的合理有限元模型，采用流固耦合方法對深圳市陰香類行道樹的破壞模式、失效風速及枝條損失率等進行系統研究，為城市防臺減災及規劃建設提供參考和指導。

1 陰香類行道樹有限元模型建立

1.1 陰香類行道樹幾何模型構建

實地調查與文獻統計顯示深圳的常見行道樹分為陰香、大葉榕、桃花心木等12種，其中陰香類所在街道的占比超過6.73%［12］，因此，本研究主要以陰香類行道樹作為研究對象進行抗風性能分析。樹木作為一種非均質材料，每個枝條的形態各異，材料參數亦各不相同，力學性能各向異性。考慮到樹木的整體與局部存在自相似、自仿射等關系，本研究將分形理論應用于陰香幾何模型的構建中，枝條在主干上分形的角度、長度和級數等參數通過分形來確定。通過對深圳大學校園內近100 多棵陰香樹的實地調研，分析統計了其主要枝條數量及生長朝向等，建模時對陰香行道樹進行了適當的簡化，本文對陰香類行道樹簡化如下：

（1）將陰香類行道樹簡化為四級體系，包含主干、一級枝條、二級枝條和三級枝條（圖1所示）；

圖1 行道樹的幾何模型Fig.1 Geometric model of border trees

（2）不考慮樹木根部與土體的耦合作用，將樹木固接在土體上；

（3）不考慮樹木初始傾斜和位移影響。

借助UG 軟件的族表和UDF 庫功能，設置Ri、Si、Ti和H0等14個系數來描述樹枝底徑、頂徑、長度、分枝位置等參數與主干底徑D0之間的關系，如圖2所示，各系數的含義見表1。考慮到耦合計算需要耗費較大的計算機資源，同時為了保證計算精度，需要控制流場域與樹木的交界面的最小網格尺寸，故計算模型根據實際樹木枝干直徑，按1：6的比例建立陰香樹的縮尺模型，其縮尺模型參數見表2。首先確定樹干和一級分枝模型，根據陰香行道樹的幾何尺寸、空間形態的特點，通過修改底徑D0和14個系數，快速獲得相應形態的樹枝模型，放置在一級分枝上，形成第二、第三級分枝，生成完整行道樹模型。

圖2 樹干及其分枝空間結構示意圖Fig.2 Spatial structure of trunk and branches

表1 模型系數及其意義Table 1 Model coefficients and their significance

表2 行道樹分形模型參數Table 2 Parameters for fractal model for border trees

樹木力學性能受其密度、樹齡及含水率等多因素影響，而且樹木不同部位力學性能亦存在差異。根據《中國木材志》可知陰香樹種與白楊樹種的基本密度和抗彎強度較為接近，樹種強度均類屬中等偏低類。因此在本研究中，陰香類行道樹的力學參數參考文獻［14］進行設置，具體參數如表3所列。另外，參考Gardinaer［13］的研究結果，一級枝條的材料參數設置為主干材料參數的80%，二級枝條為一級枝條的70%，三級枝條為二級枝條的60%，樹木的軸向彈性模量遠大于其他方向，可將行道樹模型看作橫觀各向異性體。

表3 行道樹力學性能參數Table 3 Mechanical parameters of border trees

1.2 陰香類行道樹模型固有頻率驗證

模型的準確性是風荷載及失效模式分析的基礎，國內外對針葉樹固有頻率和樹干參數的關系進行了大量的研究，Moore［15］對8個樹種602棵樹的固有頻率進行了分析，發現樹木固有頻率和胸徑（DBH）與樹高（H）的平方之比呈線性關系，采用最小二乘法擬合得如下關系式：

式中，f為樹木固有頻率；DBH為樹木的胸徑，單位cm；H為樹木的高度，單位m。

文中1.2 節所建立的行道樹簡化模型采用有限元分析得樹木的固有頻率為0.95 Hz，與采用式（1）計算得固有頻率1.09 Hz 較為接近。為了進一步驗證模型的有效性，選取5種 不 同 尺 寸（DBH/H2分 別 為0.320、0.338、0.367、0.370、0.410）的行道樹模型在ANSYS 平臺進行模態分析，以獲取模型的固有頻率。圖3所示為在全枝、去掉一個一級分枝和無枝樹干3 種情況下自振頻率f與DBH/H2的關系，圖3 中同時給出了式（1）的計算結果進行對比。由圖3 可見，行道樹固有頻率f與DBH/H2具有較好的線性關系，樹木的自振頻率隨著胸徑的增大而增大，隨著樹冠高度的增加而減小，隨著枝條的增加而減小。

圖3 自振頻率f與DBH/H2的關系圖Fig.3 Relationship between natural frequency f and DBH/H2

2 數值模擬參數設置

本研究將采用有限元分析平臺ANSYS的Fluent，Static Structural和System Coupling等模塊對行道樹的風荷載和風致振動特性進行耦合模擬分析。

2.1 計算域設置及網格劃分

本計算模型的計算域如圖4 所示，計算域長×寬×高尺寸為4 m×4 m×3 m。樹干底端位于計算域底部中心位置。計算域上風向定義為風場入口邊界、下風向定義為出口邊界，底部定義為壁面邊界、兩個側面定義為對稱邊界。樹木表面設置為光滑無滑移壁面，風場同模型外表面相接觸的面定義為流固耦合面，通過耦合面將風場產生的風荷載傳遞給行道樹模型，從而實現數據傳遞。模型最外端樹枝距計算域入口的距離為1.50 m，到出口端距離為0.96 m；距離計算域左側壁面的距離為1.33 m，右側壁面的距離為1.48 m；模型頂端距離計算域頂部的距離為0.17 m；模型疏透性較好，計算域的阻塞率為4.5%。便于簡化分析，模擬時將行道樹根端部設置為固定端約束。

圖4 計算域示意圖Fig.4 Sketch of computation domain（unit：m）

網格結構和疏密程度將直接影響到仿真模擬的計算精度，但是網格加密又會增加計算時間和儲存空間，本文在網格劃分中，合理的調整收斂控制選項，并對網格的獨立性進行了檢驗，如表4 所列。根據現有計算機能力，在保證精度的前提下，后續采用了第三套網格對行道樹進行模擬分析。

表4 25 m/s風速下的網格獨立性檢驗Table 4 Grid independence test at wind speed of 25 m/s

圖5 所示為行道樹流場網格劃分示意圖。為了提高計算精度，在流固耦合界面上設置了10層膨脹層加密網格。網格單元數量為10 010 562，扭曲度最大值為0.754，平均值為0.265；正交質量最小值為0.22，平均值為0.734，網格質量較好。為避免網格退化現象，在利用Static Structural 模塊進行受力分析時采用ANSYS 單元庫中的solid187 單元進行模擬，選擇四面體網格劃分單元，樹木實體網格單元數量為1 064 839。

圖5 網格劃分示意圖Fig.5 Sketch map of meshes

2.2 風場求解參數設置

為了在滿足收斂條件的前提下盡量節省計算資源，本文數值計算采用SIMPLE 算法對速度-壓力耦合求解；設定風場計算總時間為10 s，時間步長為0.01 s。采用標準k-ε湍流模型，每一時間步長的迭代次數設定為20 次，收斂精度設置為0.001，其他參數采用默認設置。由于縮尺模型的最大高度為2.83 m，故不考慮地面粗糙程度對風速剖面的影響，速度入口直接采用沿高度不變的均勻風速場；采用15 m/s 至35 m/s（5m/s 遞增）風速對樹木進行單向流固耦合的風荷載模擬，將模擬得到的風荷載通過流固耦合面加載到樹木模型上，進而開展樹木模型動力響應分析［16］。

3 典型行道樹風致破壞分析

3.1 行道樹風致破壞判別準則

Baker 等［17］針對云杉動力響應的研究認為樹木失效風速與結構單元平均彎矩、根端部極值彎矩以及剪切應力較為敏感，后續研究者鑒于此分別建立了基于彎矩的風致破壞準則和基于剪應力的破壞準則2類。

3.1.1 基于彎矩的破壞準則

樹木整體的風致振動主要表現為沿順風向的彎曲振動，一旦模型單元的彎矩荷載超過其抗彎能力時，即認為該單元發生彎曲破壞。樹木的最大彎矩一般出現在樹根或枝干的端部位置，因此本文以樹木根部或枝干端部的失效彎矩作為其破壞的判斷依據，當單元的最大彎矩大于其失效彎矩時，即認為發生破壞。定義彎矩的指示函數為：

式中，Ii為第i個單元的狀態指示函數，Ii＜1 意味著單元安全，Ii≥1 意味著單元失效；Ne為行道樹模型的單元數量；Bi為第i個單元的最大彎矩；B^i為第i個單元的失效彎矩。

失效彎矩B^i可根據強度理論由下式來確定［17］：

式中，σ為材料強度；ri為第i個單元的截面半徑。

文中在采用彎矩破壞準則時，首先通過強度理論計算得到枝干的失效彎矩（式（3）），之后基于風振響應分析計算出不同等級風速作用下的枝干單元截面的最大彎矩，最后由彎矩指示函數值來判定樹木破壞狀態，并對失效風速進行預測。

3.1.2 基于剪應力的破壞準則

樹木在風荷載的作用下，除了由于彎曲變形導致的彎曲破壞之外，在主干和枝條連接部位產生較大的剪切應力從而導致主干端部和枝條連接處發生剪切破壞。研究表明木材的抗剪強度僅為抗彎強度的十分之一左右［18］。因此，本文亦采用基于剪切應力的破壞準則對樹木模型的單元失效狀態進行判定［19］。

樹木的強度受生長年限、含水率、生長狀態等因素影響較大，參考鹿振友等［20］對3 年生的毛白楊的力學性能的測試結果，結合陶嗣巍對簡化模型的強度折減結果［14］，采取材料標準試驗強度值的75%作為模型的極限值，折算出樹木主干的最大剪切應力為3.75 MPa。根據Gardinaer 等［13］的研究結果，對每一級枝條的強度參數進行適當折減：一級枝條的強度參數設置為主干強度參數的80%，二級枝條為一級枝條的70%，三級枝條為二級枝條的60%。模型各分枝單元的最大切應力限值如表5所示。

表5 樹木各分枝的最大切應力Table 5 The maximum shear stress of branches of border tree

3.2 行道樹風振響應分析

3.2.1 行道樹內力響應分析

文中的模擬結果顯示樹木的最大彎矩出現在樹根端部位置，與艾曉秋等［9］的研究結果一致。圖6 所示為不同風速作用下樹根端部彎矩變化趨勢，由圖可見，隨著風速的增加，樹根端部彎矩呈冪函數遞增趨勢。

圖6 不同風速下根部彎矩Fig.6 Bending moment at the end of border tree model under different wind speeds

圖7 所示為不同風速等級下主干最大剪應力包絡圖，從圖可見，剪應力隨著平均風速的增大而增大，在最低樹枝與枝干連接處達到最大值；每個分枝受到的力，通過一級枝條傳給主干，導致每個一級枝條高度處的剪應力都會陡然增加。因此，對最低樹枝與主干的連接處采取合理的加固措施，能有效的防止樹木的風致破壞。

圖7 不同風速下主干最大切應力Fig.7 The maximum shear stress of the trunk under different wind speeds

3.2.2 行道樹位移響應分析

文中陰香類行道樹有限元模型的樹根端部采用固定支撐約束，假定其在順風向風荷載的作用下，只發生彎曲變形，不考慮其整體傾覆性破壞。在現實情況下，少量行道樹由于移植過程導致其根系較淺，在臺風過程中可能發生整體傾覆破壞。但行道樹抗傾覆力矩與樹種、根系、移植年限、土體穩定性等相關，在本研究中暫不考慮陰香行道樹的整體傾覆破壞。表6所示陰香類行道樹模型最低分枝與主干連接處在不同等級風速作用下的水平位移，以及響應的主干彎曲角度。由表可見，最低樹枝與主干連接處的位移隨著作用平均風速的增加呈逐漸增大趨勢；由最大剪切應力分析結果可知該樹木模型在風速為32 m/s 時將發生剪切破壞，此時最低樹枝與主干連接處位移約為280 mm，彎曲角度約為13.89°。

表6 行道樹模型位移Table 6 Displacements of border tree model

3.3 行道樹失效風速預測

樹木風致破壞的受力較為復雜，在考慮其受剪破壞的同時也應注意其受彎破壞的可能性。當以樹木根部彎矩作為判斷依據時，失效彎矩理論破壞系數為彎矩的指示函數；當以極限剪切應力作為失效判據時，定義破壞系數為樹干在風壓作用下樹干外層的最大應力與其相應臨界應力強度值的比。

圖8所示為行道樹根端部的彎矩破壞系數和剪力破壞系數隨風速的變化關系，由圖可見，隨著風速逐漸增大，2種破壞系數均呈逐漸增大趨勢，但是基于彎矩的破壞系數小于基于剪力的破壞系數。圖8同時給出了文獻12 的模擬結果作為對比。由圖可見，當風速達到32 m/s 時，基于最大切應力的破壞系數為100%，而基于彎矩的破壞系數則為40%。由此可見，本行道樹模型的剪切破壞是引起其結構失效的主要破壞形式，在后續分析中將主要以最大剪切應力作為主干的失效破壞準則來進行分析。

考慮到真實環境中隨著枝條的折斷掉落，樹木的迎風面積將減少，樹木所受總風荷載也將減少［21］。因此，需要針對每次模擬結果對失效枝條進行刪減，進而調整模型形態，考慮其在不同風速工況下的動力響應，以模擬枝條的動態損失過程。圖9 所示為采用最大剪應力作為枝條失效判斷依據的不同風速下樹木損失百分比。由圖可見，樹木模型當風速達到15 m/s時二級和三級枝條出現損傷掉落，考慮到三級枝條的長度和剛度較小，受到的彎矩和剪應力都較小，未發生極值失效破壞，主要隨著二級枝條的折斷而掉落，故兩者的枝條損失率相同。隨著風速增加，一級枝條隨之出現折斷掉落，導致二級和三級分枝的損失率顯著增加；同時一級枝條掉落后，導致模型整體動力特性也隨之改變，在相同風速下使得低級枝條所受風荷載降低，其破壞風速明顯提高。到達主干破壞風速（55 m/s）后，主干發生折斷破壞，樹木整體失效。與圖8對比可知，考慮枝條動態損失后，樹木模型的整體破壞風速由原來的32 m/s增加到了55 m/s，樹木的抗風性能有了顯著提高，模擬結果可為城市行道樹枝條的修剪造型提供一定的參考與借鑒意義。

圖8 破壞系數與風速關系Fig.8 The relationship between trunk damage coefficients and wind speeds

圖9 不同風速樹木損失百分比Fig.9 Percentage of tree losses under different wind speeds

4 結論

文中以深圳市陰香類行道樹作為研究對象，根據分形理論，首先建立了行道樹的簡化幾何模型；之后基于ANSYS有限元分析平臺模擬了行道樹在不同等級風速作用下的風荷載以及內力和位移，基于2種破壞準則對典型行道樹的破壞模式和失效風速進行了分析。通過本研究獲得以下主要結論：

（1）考慮行道樹的空間形態、分枝影響和力學性能差異，參數化建立了行道樹的有限元模型，該模型能合理地反映樹木的受力狀態。

（2）行道樹的彎矩和剪切應力隨著風速等級的增大而增大，其中主干的最大剪切應力出現在最低樹枝與樹干的連接處，最大彎矩則出現在樹根端部位置。

（3）樹木模型的剪切破壞是引起結構失效的主要破壞形式，彎曲破壞是次要形式。

（4）隨著枝條損失，導致模型整體動力特性也隨之改變，行道樹承受的總風荷載減少，失效風速顯著增加。