李正農(nóng)，余世斌，吳紅華，周利芬

(湖南大學(xué)建筑安全與節(jié)能教育部重點試驗室，湖南長沙，410082)

我國是風(fēng)災(zāi)頻發(fā)的國家之一，不論是沿海地區(qū)的臺風(fēng)還是西北內(nèi)陸地區(qū)的強(qiáng)風(fēng)都可能對森林生態(tài)系統(tǒng)造成災(zāi)害性影響，也會帶來巨大的經(jīng)濟(jì)損失。在強(qiáng)風(fēng)作用下，樹木可能會發(fā)生樹枝折斷甚至整體倒伏；樹木也會對風(fēng)速產(chǎn)生遮蔽和削弱效應(yīng)，改變樹木周圍流場，如降低樹后風(fēng)速、增加樹后湍流度等。研究樹木對風(fēng)場的影響對于改善近地面風(fēng)環(huán)境，減輕強(qiáng)風(fēng)對行人、建筑物或鄰近樹木的破壞等都具有十分重要的意義。目前，對于樹木周圍流場的研究，國內(nèi)外學(xué)者主要采用現(xiàn)場實測、風(fēng)洞實驗及數(shù)值模擬3種方法。在現(xiàn)場實測與風(fēng)洞實驗方面，GUAN 等[1]通過風(fēng)洞實驗，得出了氣動孔隙率、光學(xué)孔隙率以及樹木阻力系數(shù)之間的關(guān)系。CAO等[2]對不同風(fēng)速下3種喬木的阻力與傾覆力矩進(jìn)行了研究，并同時分析了樹木孔隙率、樹木迎風(fēng)面積、樹頂位移、撓度、湍流度與風(fēng)速之間的關(guān)系。MANICKATHAN 等[3]通過對比風(fēng)洞實驗中的人造模型樹與小型天然樹木之間的流場以及阻力系數(shù)的差異，發(fā)現(xiàn)只有當(dāng)模型樹和天然樹具有相似的孔隙率，且模型樹與天然樹在風(fēng)速不同情況下具備相同的空氣動力特性時，模型樹和天然樹的阻力系數(shù)才具有相似性。

在數(shù)值模擬方面，SVENSSON等[4]將樹冠簡化為高為2.5 m、阻力系數(shù)為0.3、葉面積密度為2.1 m?1的長方形來研究二維樹冠內(nèi)外的流場分布；2004—2010 年，OHASHI 等[5?7]借用附加源項法的模擬方法，將樹冠模型簡化為簡單的幾何體并等效為多孔介質(zhì)，最后以風(fēng)速為指標(biāo)深入分析了二維樹木周圍的流場。ROSENFELD 等[8]通過數(shù)值模擬并將柏樹模型簡化為三維圓錐形，以兩排(每排3 顆樹)組成三維防風(fēng)林研究了樹木之間的流場干擾作用。徐志[9]通過修正REALIZABLEk-ε模型算法，將樹冠簡化為高為5.8 m、寬為2.0 m的長方形，通過改變樹冠層疏透性、樹冠層高度、樹冠層形狀分析相應(yīng)位置風(fēng)場的變化規(guī)律。

以往的樹木周圍流場的數(shù)值模擬研究結(jié)果大多未與實測或風(fēng)洞實驗結(jié)果進(jìn)行詳細(xì)對比分析，且由于模擬技術(shù)的局限性對樹木模型進(jìn)行了一定簡化，例如簡化為長方體、圓錐體、棱臺等。本文根據(jù)本課題組所完成的樹木風(fēng)洞實驗，嚴(yán)格依據(jù)風(fēng)洞試驗時的樹木模型外輪廓建立數(shù)值模擬模型。另外，以往對樹木進(jìn)行數(shù)值模擬時，孔隙率不變，而風(fēng)洞實驗結(jié)果表明，樹木隨著風(fēng)速的增加，其孔隙率有所增大，阻力系數(shù)有所減小。此外，許多學(xué)者對樹后流場進(jìn)行分析時多只以風(fēng)速為評價指標(biāo)，對流場的湍流度分析較少。為此，在不同風(fēng)速下，本文作者采用附加源項法改變阻力系數(shù)，研究孔隙率變化對樹木的影響，并對樹木周圍湍流度進(jìn)行分析，以便為理解樹木在強(qiáng)風(fēng)作用下的破壞機(jī)理提供支持，并對樹木的種植及防護(hù)提供參考。

1 多孔介質(zhì)的數(shù)值模擬

為了研究多孔介質(zhì)的流動控制機(jī)理并進(jìn)行力學(xué)分析，需要對多孔介質(zhì)模型進(jìn)行簡化，建立合理的計算模型。現(xiàn)階段對于多孔介質(zhì)的模擬研究方法主要有2種：直接模擬法和附加源項法。直接模擬法是直接按照多孔介質(zhì)對樹木模型進(jìn)行一定數(shù)量的開孔，再對其進(jìn)行網(wǎng)格劃分，這種方法因模型過于復(fù)雜較難實現(xiàn)，開孔區(qū)域模型的復(fù)雜性導(dǎo)致網(wǎng)格在該區(qū)域易發(fā)生畸變，從而影響數(shù)值計算效率與精度。附加源項法是基于實驗測得氣流流經(jīng)多孔介質(zhì)區(qū)域后速度、湍流度等一系列物理參數(shù)的變化情況，通過編譯UDF(user define function)的形式，對多孔介質(zhì)流域的動量、湍動能以及湍動能耗散率方程進(jìn)行修正，通過在動量守恒方程、湍動能和湍動能耗散率方程中添加附加源項，對氣流流經(jīng)樹木后的流動狀況進(jìn)行模擬。

對于厚度不是很大的多孔介質(zhì)如很薄的擋風(fēng)板，氣流流過該區(qū)域后一般不會發(fā)生很復(fù)雜的流動，因而只對其動量方程添加附加源項即可。動量方程修正后的表達(dá)方式如下：

式中：Si為量方程的附加源項，定義為氣流流經(jīng)多孔介質(zhì)域損失量；Dij為阻力系數(shù)矩陣；Cij為阻力系數(shù)矩陣；μ為氣流的黏性系數(shù)；vi為氣流在i方向上的速度；ρ為氣流的密度，取1.225 kg/m3。

由于實際上氣流雷諾數(shù)偏大，相較于其慣性阻力系數(shù)，其黏性阻力系數(shù)可以忽略不計。同時，為了簡化計算，只需考慮氣流在流動方向上的速度，Si簡化為

式中：Cf為單位長度的阻力系數(shù)；Cd為樹冠層阻力系數(shù)；η為樹冠層葉面積指數(shù)(一般樹種取2～4 m?1)，即單位土地面積上植物葉片總面積占土地面積的倍數(shù)[9]；D為多孔介質(zhì)區(qū)域厚度，m；F為樹冠層所受的阻力，N；Uref為參考風(fēng)速，m/s；A為樹木靜止時迎風(fēng)面的投影面積，m2。

樹冠對湍流的影響由標(biāo)準(zhǔn)的k-ε模型修改。修改的基礎(chǔ)是在樹冠層內(nèi)部的湍動能方程和湍動能耗散率方程中加入源項，可進(jìn)行如下表達(dá)：

式中：μτ為湍流黏性系數(shù)；Gk為平均速度梯度引起的湍動能；常數(shù)C1ε=1.44，C2ε=1.92；σk為湍動能的湍流普朗特數(shù)；Sk為湍動能附加源項表達(dá)式；Sε為湍動能耗散率附加源項表達(dá)式。對于模擬近地面氣流流經(jīng)樹木后的湍動能和湍動能耗散率附加源項模型Sk和Sε，國外學(xué)者進(jìn)行了相關(guān)研究，國內(nèi)對這方面的研究較少。數(shù)值模型大概可分為4種類型，見表1。

根據(jù)MOCHIDA等[10]簡化的樹木模型，并在前人研究的基礎(chǔ)上，發(fā)現(xiàn)上述模型的第4 組(見表1)更能有效地模擬樹木對湍流的影響，可將Sk和Sε分別表示為：

表1 湍動能和湍動能耗散率附加源項Table 1 Additional source terms for turbulent kinetic energy and turbulent kinetic energy dissipation rate

式中：Cρε1和Cρε2為附加源項系數(shù)，用于控制多孔介質(zhì)中能量耗散過程的時間尺度。根據(jù)MOCHIDA等[10]算例中附加源項法各項系數(shù)的取值，建議將Cρε1和Cρε2分別取為1.8和1.5。

2 計算模型與數(shù)值方法

2.1 計算模型與網(wǎng)格劃分

本文在計算時，樹木模型取自風(fēng)洞試驗所用的豆瓣黃楊幼苗，該樹木具有與實際樹木相似的樹冠形態(tài)、分支結(jié)構(gòu)的材料構(gòu)成，且樹干與樹枝、樹枝與樹葉之間的連接構(gòu)造與真實樹木的相同，因而，可以較好地模擬樹木在風(fēng)場中的反映。具體樹木模型見圖1。

圖1 風(fēng)洞實驗樹木模型Fig.1 Wind tunnel experiments with tree models

本文的樹木模型高度H為50.0 cm，樹干高度為10.0 cm，樹冠最大厚度為25.4 cm，各截面半徑見表2。

表2 樹木各截面半徑Table 2 Radius of each section of tree

該樹木由于近似呈現(xiàn)中心對稱圖形，因此，在SolidWorks 利用相關(guān)點坐標(biāo)，將1/2 縱截面通過旋轉(zhuǎn)360°形成3D模型，見圖2。

圖2 數(shù)值風(fēng)洞樹木模型Fig.2 Model of numerical wind tunnel tree

計算域長×寬×高為750 cm×500 cm×250 cm，阻塞率為0.7%(<3%)，采用Icem 軟件生成樹木整體網(wǎng)格。由于樹木形狀不規(guī)則，因此，采用內(nèi)外兩部分，內(nèi)部帶有樹木的網(wǎng)格采用非結(jié)構(gòu)網(wǎng)格，其區(qū)域長×寬×高為20 cm×20 cm×60 cm；外部網(wǎng)格采用結(jié)構(gòu)網(wǎng)格；總網(wǎng)格數(shù)目為202.7萬個。具體網(wǎng)格分布見圖3。

圖3 網(wǎng)格分布Fig.3 Distribution of mesh

2.2 邊界條件

計算域及邊界條件見圖4。計算域左邊界設(shè)置成入口(inlet)，右邊設(shè)置成自由出口(outflow)，外網(wǎng)格四周設(shè)置為無滑移的靜止壁面(wall)，內(nèi)外網(wǎng)格交界面設(shè)置為兩流域交界面(interface)，樹冠設(shè)置為內(nèi)部面(interior)。

圖4 計算域及邊界條件Fig.4 Computational domain and boundary conditions

入口與風(fēng)洞實驗一致的B類地貌下的風(fēng)場，其粗糙度系數(shù)α=0.15，縮尺比為1∶40，風(fēng)剖面采用指數(shù)率分布，具體風(fēng)剖面公式以及湍動能、湍動能耗散率公式見式(10)，(11)和(12)。風(fēng)場的平均風(fēng)速剖面和湍流度剖面見圖5和圖6。從圖5和圖6可以看出：模擬風(fēng)場與實際風(fēng)場較吻合，并且模擬湍流度與實際湍流度也較匹配，這在一定程度上驗證了風(fēng)場的穩(wěn)定性。

圖5 左邊界速度進(jìn)口模擬風(fēng)速與實測風(fēng)速對比Fig.5 Comparison of simulated wind speed and measured wind speed at the inlet of left boundary speed

圖6 左邊界速度進(jìn)口模擬湍流度與實測湍流度對比Fig.6 Comparison of simulated turbulivity and measured turbulivity at the inlet of left boundary speed

式中：z為高度；z0為標(biāo)準(zhǔn)參考高度；U0為標(biāo)準(zhǔn)參考高度下的平均風(fēng)速；U為該高度處的平均風(fēng)速；k為湍動能；Iu為湍流強(qiáng)度；Cμ為常數(shù)；ε為湍動能耗散率；l為湍流積分尺度。

根據(jù)何定穎[11]的實驗，分別以控制風(fēng)速8，10，12和14 m/s，表示不同風(fēng)速下的B類風(fēng)場。具體控制風(fēng)速對應(yīng)實際離地高度10 m 處(風(fēng)洞25cm高度處)的風(fēng)速見表3。

表3 控制風(fēng)速Ucon與樹冠風(fēng)速U10對應(yīng)關(guān)系Table 3 Control wind speed and corresponding canopy wind speed

2.3 網(wǎng)格無關(guān)性的驗證

現(xiàn)有模型的所有要素如冠層的多孔模型、湍流模型(包括修正的源項)和計算模型的范圍都在以往的研究中得到了檢驗和驗證[12?14]，特別是Fluent CFD 軟件包的準(zhǔn)確性和適用性得到廣泛測試。由于本研究中使用的幾何形狀為樹木的外輪廓，為了驗證模型網(wǎng)格，本文對網(wǎng)格的獨立性進(jìn)行研究。

為了對比網(wǎng)格數(shù)目對計算結(jié)果的影響，采用202.7 萬、411.5 萬和963.1 萬個網(wǎng)格(見圖7)進(jìn)行計算。網(wǎng)格總數(shù)調(diào)整通過改變整體網(wǎng)格數(shù)目實現(xiàn)。

圖7 順風(fēng)向網(wǎng)格布置Fig.7 Distribution of chord mesh

由實測風(fēng)洞實驗控制風(fēng)速為10 m/s，阻力系數(shù)Cd為0.794，利用式(4)求得單位長度的阻力系數(shù)Cf為10.8，利用UDF導(dǎo)入Fluent中得到不同網(wǎng)格在高度1/2H處的軸向速度分布，見圖8。

由圖8可知：為了使數(shù)據(jù)變化更明顯，以高度0.25 m(即1/2H)處沿順風(fēng)向取點，該點處于樹冠影響最大的位置。雖然網(wǎng)格呈倍數(shù)增大，但對應(yīng)求得的數(shù)值解并沒有較大差別，最大相對誤差為0.8%，說明模擬結(jié)果對網(wǎng)格數(shù)量的靈敏度較低，數(shù)值解較準(zhǔn)確。

圖8 不同網(wǎng)格在高度1/2H軸向速度分布Fig.8 Axial velocity profiles of different grids at height of 1/2H

2.4 流場計算方法

采用Fluent作為計算工具，求解流動控制方程獲得樹木流場的數(shù)值解。

通過穩(wěn)態(tài)計算求解平均風(fēng)速時，湍流模型選用具有比標(biāo)準(zhǔn)的k-ε模型精度高和穩(wěn)定性強(qiáng)的RNGk-ε模型，并使用標(biāo)準(zhǔn)壁面函數(shù)(standard wall functions)對近壁面進(jìn)行處理。流體材料選用理想不可壓縮氣體，材料參數(shù)使用默認(rèn)值。采用SIMPLEC 方法對速度壓力方程組進(jìn)行解耦，SIMPLIC算法與SIMPLE算法相比，計算的收斂速度加快。離散格式使用二階迎風(fēng)格式(second order upwind)。穩(wěn)態(tài)模擬計算的收斂標(biāo)準(zhǔn)為各項殘差下降至10?4以下。

通過瞬態(tài)計算求解脈動風(fēng)速，湍流模型采用可精確求解湍流尺度的大渦模擬(LES)。在穩(wěn)態(tài)計算收斂后，以穩(wěn)態(tài)計算結(jié)果作為初始流場，從而加快大渦模擬的收斂速度。進(jìn)行大渦模擬計算時，使用WALE 模型，時間步長設(shè)置為Δt=0.005 s，時間步數(shù)為4 000步，子步迭代數(shù)為20。

冬季和凍融期生物過程對N2O排放具有重要影響的報道最早出現(xiàn)在1998年[36]。在寒冷冬季，深層土壤仍存在微生物活動，且積雪覆蓋維持了雪下土壤相對溫暖和穩(wěn)定的環(huán)境，導(dǎo)致冬季土壤呼吸及CO2的持續(xù)排放[20，37]。Fahnestock et al.發(fā)現(xiàn)亞北極苔原在非生長季有土壤CO2排放[38]，并由此可以估算出亞北極苔原冬季CO2總排放量在全年占一定地位。更進(jìn)一步分析，如果將這些地區(qū)冬季溫室氣體的排放量納入其全年的C收支，將使亞苔原生態(tài)系統(tǒng)、高緯度地區(qū)生態(tài)系統(tǒng)的年C通量有不同程度的增加，甚至使其由溫室氣體的“匯”轉(zhuǎn)變?yōu)椤霸础盵39，40]。

3 計算結(jié)果與分析

3.1 10 m/s 控制風(fēng)速下樹后風(fēng)剖面的模擬結(jié)果與實驗結(jié)果對比

首先根據(jù)模型的樹冠參數(shù)，將樹冠直徑10～50 cm按照所取坐標(biāo)點分成10部分，利用面積法求解樹冠的平均厚度；根據(jù)實驗結(jié)果，利用式(4)得知10 m/s 控制風(fēng)速下黃楊樹苗的試驗參數(shù)，見表4。

表4 黃楊數(shù)值模擬的參數(shù)Table 4 Numerical parameters of Boxwood

利用上述附加源項法將多孔介質(zhì)相關(guān)系數(shù)添加入Fluent 的UDF 中。為了更好地分析模擬結(jié)果與實驗結(jié)果的擬合程度，將樹后不同距離的剖面與實驗結(jié)果進(jìn)行對比，見圖9和圖10。

圖9 樹后不同距離剖面相對參考點風(fēng)速νi/νref模擬結(jié)果與實驗結(jié)果對比Fig.9 Simulation and experimental comparison of relative reference point wind speed in different distance profiles behind trees

為了方便比較，將模擬結(jié)果進(jìn)行歸一化處理，即將各測點的平均風(fēng)速除以參考點處的平均風(fēng)速，各測點的湍流度除以參考點處的湍流度，即

其中：vi為測點i處的平均風(fēng)速；vref為參考點高度處的平均風(fēng)速；Ii為測點i處的湍流度；Iref為參考點高度處的湍流度；σ為順風(fēng)向脈動風(fēng)速均方根；uˉ為順風(fēng)向平均風(fēng)速。

根據(jù)圖像觀察，0.2H～H是樹后風(fēng)場降低幅度最大的區(qū)間，這也是樹冠所在的區(qū)間，與實驗結(jié)果一致，說明實驗與模擬的樹后整體趨勢一致。模擬中相對參考點風(fēng)速最小值在最靠近樹木的X=0.6H處的樹冠中心位置，其值為0.223，而實驗中相對參考點風(fēng)速最小值為0.245，同樣在樹冠中心的位置，說明模擬結(jié)果與實測結(jié)果在數(shù)值上差別不大，而且在2H處，樹后風(fēng)速恢復(fù)了60%，這與實驗結(jié)果相符。

從圖9 和圖10 可知，0.2H～H是樹后湍流度影響幅度最大的區(qū)間，這也是樹冠所在的區(qū)間，與實驗結(jié)果一致，說明實驗與模擬的樹后風(fēng)場整體趨勢一致。模擬中相對參考點湍流度最大值在最靠近樹木的X=0.6H處的樹冠中心位置，其值為1.76，而實驗中相對參考點湍流度最大值為1.65，同樣在樹冠中心位置，說明模擬結(jié)果與實測結(jié)果差別不大，而且在2H處，樹后湍流度恢復(fù)到樹前湍流度的80%，這與實驗結(jié)果相符。

圖10 樹后不同距離剖面相對參考點湍流度Ii/Iref模擬結(jié)果與實驗結(jié)果對比Fig.10 Simulation and experimental comparison of turbulence degree relative to reference point in different distance profiles behind trees

3.2 10 m/s控制風(fēng)速下樹后特定平面流場分析

為便于分析，以下結(jié)果均采用相對風(fēng)速和相對湍流度進(jìn)行分析，即將各測點處的平均風(fēng)速和湍流度值除以無樹時流場中相對應(yīng)高度處平均風(fēng)速νrel和湍流度Irel，計算公式如下：

其中：vi為測點i處平均風(fēng)速；v0i為無樹時流場中測點i處的平均風(fēng)速；Ii為測點i處湍流度；I0i為無樹時流場中測點i處湍流度。

3.2.1 橫風(fēng)向平面相對風(fēng)速和相對湍流度分析

X=H時單棵樹相對風(fēng)速等值線見圖11。從圖11 可見：在X=H的橫風(fēng)向平面內(nèi)，樹木的相對速度以樹冠中心(Y/D=0，Z/H=0.6)呈放射狀向四周逐次增大，在樹冠中心達(dá)到最小值0.5，數(shù)值結(jié)果與風(fēng)洞實驗結(jié)果一致；在水平方向上，當(dāng)?0.5DH(樹高)時，相對速度也大于0.9，因而可以忽略樹木對于風(fēng)速的影響；當(dāng)Z<0.2H時，由于樹葉在樹冠位置，而在0.2H之下能影響流場的只有樹干，因而，相比于0.2H～H處稠密度大、葉面積密度大的樹冠區(qū)域影響要小得多。以上結(jié)果與風(fēng)洞試驗結(jié)果相符。

圖11 X=H時單棵樹相對風(fēng)速等值線圖Fig.11 Relative wind speed contour map of a single tree when X=H

用相對來流風(fēng)速減少20%即相對風(fēng)速為0.8所對應(yīng)的距離d20描述樹木的遮蔽效應(yīng)[18]。Y/D=±0.5為樹木的外輪廓范圍，可發(fā)現(xiàn)d20的形狀大致與樹冠的輪廓線一致，但在靠近地面的區(qū)域，由于樹干的寬度遠(yuǎn)小于樹冠的寬度，對風(fēng)速的減弱效果不明顯，即樹冠的形狀決定了樹木樹后橫風(fēng)向平面的有效遮蔽范圍，在風(fēng)洞實驗中也得出了類似結(jié)論。

X=H處單棵樹相對湍流度等值線見圖12。從圖12 可以發(fā)現(xiàn)：在X=H的橫風(fēng)向平面內(nèi)，樹木的相對湍流度以樹冠中心(Y/D=0，Z/H=0.6)呈放射狀向四周逐漸減小，在樹冠中心達(dá)到最大值2.49，且等值線的形狀與樹冠外輪廓線相似，這與風(fēng)洞實驗結(jié)果相一致；在水平方向上，當(dāng)距離達(dá)到±0.75D時，相對湍流度達(dá)到1.42 且在這個范圍之外，相對湍流度是逐漸接近于1；在高度方向，當(dāng)Z>H(樹高)時，向上則相對湍流度接近于1，可以認(rèn)為在樹高范圍外樹木對樹后流場的影響可以忽略不計；當(dāng)Z<0.2H時，向下則相對湍流度逐漸減小，在0.2H位置，由于其下僅有樹干存在，對風(fēng)場湍流度影響較小，相對湍流度接近于1.0，模擬中樹木對湍流度的影響與風(fēng)洞實驗結(jié)果相一致。

圖12 X=H時單棵樹相對湍流度等值線圖Fig.12 Relative turbulence contour map of a single tree when X=H

3.2.2 順風(fēng)面、水平面相對速度和湍流度分析

單棵樹水平剖面風(fēng)速云圖見圖13。由圖13 可知：在水平剖面，風(fēng)速減弱范圍為橢圓形區(qū)域，越靠近樹木，樹木對風(fēng)速的削弱作用越大，且在樹冠中心相對風(fēng)速達(dá)到最小值。結(jié)合圖14，在垂直剖面中，3H處樹后速度恢復(fù)樹前速度80%，相當(dāng)于樹木的遮蔽區(qū)域d20在X=3H以內(nèi)；此外，隨著順風(fēng)向的距離增加，相鄰兩等值線之間的距離也逐漸增大，即在距離樹木越遠(yuǎn)之處，降低同樣的風(fēng)速所需要的距離也越長；在高度方向上，可以明顯地觀察到在樹冠高度以內(nèi)，風(fēng)速減弱效果明顯，在樹冠高度以外，相對風(fēng)速接近于1，即對風(fēng)速的影響較小。

圖13 單棵樹水平剖面風(fēng)速云圖Fig.13 Nephogram of single tree horizontal profile wind speed

圖14 樹后Y=0平面相對風(fēng)速等值線圖Fig.14 Relative wind speed contour map behind tree in Y=0 plane

樹后Y=0 平面相對湍流度等值線見圖15。由圖15 可知：在Y=0 的順風(fēng)向平面內(nèi)，樹木的相對湍流度以樹冠中心呈半橢圓線向外逐漸減小，且在X=0.5H的樹冠為中心處相對湍流度達(dá)到最大值2.67，越靠近樹木，樹木對湍流度的影響越大，在樹冠中心相對湍流度達(dá)到最小值，這與風(fēng)洞實驗結(jié)果相同；在2H處，樹后湍流度恢復(fù)樹前湍流度的80%，且隨著順風(fēng)向的距離增加，相鄰兩等值線之間的距離也逐漸增大，即距離樹木越遠(yuǎn)，降低同樣的湍流度所需要的距離越長；在高度方向上，可以明顯地觀察到在樹冠高度以內(nèi)，湍流度增強(qiáng)效果明顯。

圖15 樹后Y=0平面相對湍流度等值線圖Fig.15 Relative turbulence contour map in the plane Y=0 behind tree

3.3 不同控制風(fēng)速下樹后剖面分析

該風(fēng)洞實驗通過測量風(fēng)速為8，10，12 和14 m/s時的流場，研究風(fēng)速對樹木形態(tài)的影響。利用何定穎等[12]實驗結(jié)果，得到光學(xué)孔隙率β以及實驗所測的阻力系數(shù)見表5，光學(xué)孔隙率β與單位長度的阻力系數(shù)Cf之間的關(guān)系見圖16。

表5 實測光學(xué)孔隙率和阻力系數(shù)Table 5 Measured optical porosity and resistance coefficient

圖16 光學(xué)孔隙率與阻力系數(shù)的關(guān)系Fig.16 Relationship between optical porosity and resistance coefficient

擬合的公式以及相關(guān)系數(shù)為

采用指冪指數(shù)的形式擬合光學(xué)孔隙率與阻力系數(shù)之間的關(guān)系，發(fā)現(xiàn)擬合的相關(guān)系數(shù)R2接近于1，說明擬合優(yōu)度較高，再根據(jù)表2 以及圖16，發(fā)現(xiàn)在同一樹種下，隨著風(fēng)速增大，樹木孔隙率逐漸增大，同時，阻力系數(shù)隨著孔隙率增大而減小。

在不同控制風(fēng)速下，利用不同阻力系數(shù)得到對應(yīng)的相對參考點風(fēng)速剖面以及相對參考點湍流度剖面，分別見圖17和圖18。

圖17 不同風(fēng)速下相對參考點風(fēng)速剖面Fig.17 Relative reference point wind speed profile at different wind speeds

由圖17 可知：在各個控制風(fēng)速下，數(shù)值模擬和風(fēng)洞實驗得到的樹后相對參考點風(fēng)速剖面趨勢基本一致，即在相對高度為(0～0.6)H時，隨著高度增加，相對參考點風(fēng)速呈減小趨勢，在相對高度為(0.6～1.4)H時，隨著高度增加，相對參考點風(fēng)速呈增大趨勢，并且在0.6H左右達(dá)到最小值。0.6H為樹冠中心對應(yīng)的高度，說明樹冠中心對風(fēng)速的減弱效果最明顯；此外，在樹冠高度范圍內(nèi)，隨著控制風(fēng)速增大，相對參考點風(fēng)速剖面逐漸向右偏移，相對參考點風(fēng)速最大值由8 m/s 時的0.39 增大到14 m/s時的0.47。這是因為隨著風(fēng)速增大，樹冠較靈活的樹葉和柔軟的樹枝受風(fēng)速影響，沿著風(fēng)向樹葉和樹枝進(jìn)行了重新定向，樹冠的變形會影響流經(jīng)樹木的氣流流量以及由此產(chǎn)生的阻力系數(shù)。這種現(xiàn)象在成熟的樹木現(xiàn)場測量中也有類似發(fā)現(xiàn)[15?17]。

不同風(fēng)速下相對參考點湍流度剖面見圖18。由圖18 可知：在相對高度處于(0～0.6)H時，隨著高度增加，相對參考點湍流度呈現(xiàn)增大趨勢，在相對高度為(0.6～1.4)H時，隨著高度增加，相對參考點湍流度呈減小趨勢，并且在樹冠中心處(0.6H)達(dá)到最大值；與風(fēng)速分析結(jié)果相似，隨著控制風(fēng)速增大，樹冠收縮變形，氣流在通過樹冠時受到樹枝和樹葉的擾動相對減弱，體現(xiàn)在樹冠高度范圍內(nèi)相對參考點隨著風(fēng)速的增大而逐漸向左偏移。相對參考點湍流度的最大值由8 m/s 時的3.53 減小到14 m/s 時的2.79。同理，也是因為隨著風(fēng)速增加，樹木樹冠形狀發(fā)生了改變，使得孔隙率增大，對湍流度的影響也逐漸減小。

圖18 不同風(fēng)速下相對參考點湍流度剖面Fig.18 Relative reference point turbulence profile at different wind speeds

4 結(jié)論

1)利用Fluent可以有效模擬樹木對風(fēng)場的干擾效果，而且可以通過定量阻力系數(shù)模擬不同風(fēng)速對樹木周圍風(fēng)場的影響效果。

2)樹木對樹后流場的影響主要以樹冠為中心，以橢圓形呈輻射狀想四周擴(kuò)展，距離樹冠中心越遠(yuǎn)，樹木對風(fēng)速的削弱效果越小，對湍流度的影響也越小。在樹后橫風(fēng)向平面(X=H)內(nèi)，樹木的遮蔽區(qū)域與樹冠的區(qū)域相重合，與風(fēng)洞試驗結(jié)果一致。

3)對于同一種樹木，隨著風(fēng)速增加，樹木孔隙率逐漸增大，阻力系數(shù)相應(yīng)減?。痪嚯x樹木1倍樹高處，風(fēng)剖面的相對參考點風(fēng)速最小值由0.39增大到0.47，相對參考點湍流度最大值由3.53減小至2.79；考慮因風(fēng)速增大而導(dǎo)致樹木阻力系數(shù)減小的情況，模擬的風(fēng)場與風(fēng)洞試驗測量的風(fēng)場仍然類似，說明附加源項法可以較好地模擬樹后風(fēng)場。