樹冠流動(dòng)阻力特性數(shù)值模擬與實(shí)驗(yàn)研究

楊會(huì)，付海明

?

樹冠流動(dòng)阻力特性數(shù)值模擬與實(shí)驗(yàn)研究

楊會(huì)1, 2，付海明1

(1. 東華大學(xué)環(huán)境科學(xué)與工程學(xué)院，上海，200051；2. 桂林航天工業(yè)學(xué)院能源與建筑環(huán)境學(xué)院，廣西桂林，541004)

為深入分析樹冠內(nèi)部結(jié)構(gòu)參數(shù)對(duì)樹冠繞流流動(dòng)阻力的影響，采用數(shù)值模擬和實(shí)驗(yàn)測試相結(jié)合的方法研究幾種典型形態(tài)特征樹冠的內(nèi)部流場特征和流動(dòng)阻力。用樹冠葉面積指數(shù)(AI)、稠密度(VF)和分形維數(shù)(f)表征樹冠復(fù)雜的形態(tài)結(jié)構(gòu)，討論樹冠形態(tài)結(jié)構(gòu)參數(shù)對(duì)樹冠流動(dòng)阻力的影響規(guī)律。研究結(jié)果表明：模型葉面積指數(shù)與稠密度呈近似線性關(guān)系，而樹冠分形維數(shù)分別與樹冠稠密度和葉面積指數(shù)滿足近似對(duì)數(shù)關(guān)系。通過對(duì)數(shù)值模擬計(jì)算結(jié)果回歸分析，給出樹冠流動(dòng)阻力系數(shù)與葉面積指數(shù)的關(guān)系式。研究結(jié)果表明：真實(shí)條件下的樹冠流動(dòng)阻力結(jié)果要明顯高于二維簡化模型的數(shù)值計(jì)算結(jié)果，但二者在變化趨勢上表現(xiàn)出較好的一致性；在引入樹冠阻力系數(shù)修正系數(shù)后，數(shù)值計(jì)算結(jié)果與實(shí)驗(yàn)結(jié)果吻合較好。

樹冠；流動(dòng)阻力；樹冠阻力系數(shù)；葉面積指數(shù)

樹具有保持水土、防風(fēng)固沙、涵養(yǎng)水源、調(diào)節(jié)氣候、減少污染物等功能。樹冠系指喬木樹干和其他的連同集生枝葉的部分，其狀如冠，通常將自然形成的樹冠稱為自然樹冠，而將經(jīng)過人工修剪造型而成的樹冠則稱為人工樹冠。長期以來，針對(duì)風(fēng)場作用下樹冠流場特性研究主要是采用野外觀測、風(fēng)洞試驗(yàn)以及數(shù)值模擬。如CAO等[1]通過邊界層風(fēng)動(dòng)實(shí)驗(yàn)裝置測定3種樹冠阻力特性、風(fēng)速分布以及尖端位移; SANZ等[2]于風(fēng)洞實(shí)驗(yàn)中使用粒子成像測速儀測定二維森林冠層流動(dòng)及湍流特性；MOLINA-AIZ等[3]通過風(fēng)洞試驗(yàn)測了4種園藝植物的阻力系數(shù)D；DUMAN等[4]在大氣邊界層內(nèi)實(shí)測得了樹冠氣流的速度損耗和湍流耗散；IRVINE等[5]則實(shí)測了北美典型云杉林的風(fēng)場和湍流強(qiáng)度特征。NOVAK等[6]利用Dantec三軸熱線探針測定變?nèi)~面積密度的人造樹的風(fēng)速及湍動(dòng)能分布。國內(nèi)也有少量研究者開展了樹冠流場特性的研究，如關(guān)德新等[7]利用多年的林帶疏透度、透風(fēng)系數(shù)的野外觀測和風(fēng)洞實(shí)驗(yàn)測試結(jié)果給出了立體型林帶和平面型林帶疏透度和透風(fēng)系數(shù)的關(guān)系式，在此基礎(chǔ)上討論了林帶阻力系數(shù)與疏透度和透風(fēng)系數(shù)的關(guān)系。目前，由于測試條件的限制，試圖通過實(shí)驗(yàn)手段完全了解樹冠內(nèi)外流場分布仍十分困難。如野外對(duì)真實(shí)樹冠進(jìn)行流場測試時(shí)，來流風(fēng)速不受人為控制，且樹冠較為高大不易測量；若進(jìn)行風(fēng)洞測試，則因風(fēng)洞較小，難以獲得真實(shí)環(huán)境下近地面來流風(fēng)場速度分布；此外，實(shí)驗(yàn)也難以獲得樹冠內(nèi)外流場的精細(xì)結(jié)構(gòu)。數(shù)值模擬則提供了一種有效的研究方法。ZENG等[8]利用一階閉合模型來預(yù)測樹冠內(nèi)外空氣流速，而SIQUEIRA等[9]的研究結(jié)果表明二階閉合模型在求解輸運(yùn)方程方面更有優(yōu)勢。SOGACHEV等[10?11]提出并驗(yàn)證考慮樹冠阻力的改進(jìn)兩方程(?)模型以及一種通過增加附加方程來減小模擬擬冠層流動(dòng)和邊界層流動(dòng)不確定的方法。MOCHIDA等[12]將源于實(shí)驗(yàn)結(jié)果的模型參數(shù)運(yùn)用?湍流模型，并預(yù)測風(fēng)速分布以及能量損失。HIRAOKA等[13]在假設(shè)彌散量遠(yuǎn)小于雷諾應(yīng)力的前提下運(yùn)用?湍流模型模擬樹冠流動(dòng)，所得數(shù)值模擬結(jié)果與實(shí)驗(yàn)結(jié)果吻合良好。李亮等[14]用帶源項(xiàng)模型模擬了三維樹冠流場分布，結(jié)果表明數(shù)值計(jì)算獲得的速度、及參量的空間分布與GREEN[15]的實(shí)驗(yàn)值吻合較好。HUANG等[16?17]運(yùn)用大渦模型(LES)模擬三維森林邊界湍流特性。DUPONT等[18]運(yùn)用大渦模型(LES)以及改進(jìn)的氣象模型高級(jí)區(qū)域預(yù)報(bào)系統(tǒng)模擬樹冠內(nèi)外湍流特性。PATTON等[19]利用大渦模型(LES)模擬并分析覆蓋有稀疏及濃密斜坡的壓力場和速度分布。YANG 等[20?21]也采用大渦模型分析了空氣穿越樹冠時(shí)所產(chǎn)生的回流特征，并討論了回流區(qū)與葉面積指數(shù)(AI)之間的相互關(guān)系；BELCHER等[22]提出綜合分析討論復(fù)雜地勢和森林邊界對(duì)湍流統(tǒng)計(jì)量的影響。然而，上述研究主要是針對(duì)樹冠外部的流場特征，且未涉及樹冠阻力系數(shù)與樹冠形態(tài)特征的關(guān)系。由于樹冠的樹葉形式的多樣性及樹葉數(shù)量巨大，直接模擬具有葉片尺度特征的樹冠內(nèi)部流場的精細(xì)結(jié)構(gòu)目前尚存諸多計(jì)算手段上的困難，不同形狀的樹冠對(duì)流場的影響也有所不同，對(duì)樹冠所作簡化越多，相應(yīng)的計(jì)算結(jié)果誤差也會(huì)偏大。SVENSSON等[23]將樹冠簡化為2.5 m高、阻力系數(shù)為0.3、葉面積密度為2.1 m?1的長方形來研究二維樹冠內(nèi)外的流場分布；同樣地，IRVINE等[5]和GASH[24]也分別把三維樹冠簡化為二維長方形區(qū)域來研究森林邊界處流場分布。李亮等[14]將三維樹冠簡化為四棱錐、四棱臺(tái)、圓錐及長方體，比較分析不同簡化的樹冠模型內(nèi)外流場分布差異。在前人的基礎(chǔ)上，本文作者考慮將三維復(fù)雜的樹冠結(jié)構(gòu)進(jìn)行簡化，采用圓環(huán)表征由眾多樹葉組成的團(tuán)簇，建立二維樹冠結(jié)構(gòu)簡化模型。在樹葉團(tuán)簇尺度上模擬樹冠內(nèi)部流場的精細(xì)結(jié)構(gòu)和阻力特性。為獲得真實(shí)條件下樹冠的流動(dòng)阻力特性，選取4種典型形態(tài)結(jié)構(gòu)的樹枝進(jìn)行風(fēng)洞試驗(yàn)，進(jìn)一步探索樹冠阻力特性與其形態(tài)結(jié)構(gòu)的關(guān)聯(lián)性，試圖給出樹冠流動(dòng)阻力系數(shù)估計(jì)的簡潔模型。

1 物理模型

考慮圖1所示的二維計(jì)算平面域，計(jì)算流域高為9 m，流域長為39 m。樹冠高度為3 m，流域高度設(shè)為9 m，即三倍樹冠高度，故可忽略流域高度對(duì)流場求解結(jié)果的影響。為保證入口邊界和出口邊界處流動(dòng)的充分發(fā)展，設(shè)定樹冠到流體入口距離為9 m，亦為3倍樹高，樹冠到流體出口邊界距離為30 m，即10倍樹高。流域被劃分為390×90個(gè)非均勻三角形網(wǎng)格，并對(duì)樹冠內(nèi)部作局部網(wǎng)格加密處理。邊界條件設(shè)置如下：地面設(shè)置為無滑移粗糙壁面，粗糙度取5 mm，并假設(shè)壁面無明顯的障礙物和植被；流域上部設(shè)置為滑移邊界，滑移速度與來流速度相同；流體進(jìn)出口分別設(shè)為速度入口邊界和壓力出口邊界；樹冠輪廓和樹葉團(tuán)簇輪廓分別設(shè)置為內(nèi)部邊界(interior邊界)和粗糙無滑移壁面，粗糙度取0.2 mm。

為簡化分析，現(xiàn)對(duì)問題進(jìn)行如下基本假設(shè)：1) 流體流動(dòng)為不可壓縮、無旋、定常二維流；2) 大氣環(huán)境為中性，且忽略Coriolis力，其中Coriolis力是對(duì)旋轉(zhuǎn)體系中進(jìn)行直線運(yùn)動(dòng)的質(zhì)點(diǎn)由于慣性相對(duì)于旋轉(zhuǎn)體系產(chǎn)生的直線運(yùn)動(dòng)的偏移的一種描述，源于物體運(yùn)動(dòng)所具有的慣性。

圖1 計(jì)算區(qū)域及邊界條件(樹葉團(tuán)簇的當(dāng)量直徑DI=4.84 cm，LAI=5.3)

在上述假設(shè)條件下，樹冠繞流流動(dòng)的動(dòng)量方程和連續(xù)性方程分別為：

式中：U為方向上平均速度分量；為實(shí)際壓力與參考值之差；u為平均速度U的脈動(dòng)值；0為空氣密度；為空氣運(yùn)動(dòng)學(xué)黏度。

式中：μ為渦黏性；δ為Kronecker函數(shù)；為湍流動(dòng)能。由和定義速度的比例尺和長度比例尺分別為

在湍流模型中，渦黏性μ定義為

式中：C為經(jīng)驗(yàn)常數(shù)。

本文采用LAUNDER等[25]提出的標(biāo)準(zhǔn)?湍流模型。湍動(dòng)能及其耗散速率方程分別為：

式中：σ，σ，C1和C2為經(jīng)驗(yàn)常數(shù)[26]，見表1。

表1 標(biāo)準(zhǔn)k?ε湍流模型各經(jīng)驗(yàn)常數(shù)

2 風(fēng)洞實(shí)驗(yàn)

為了驗(yàn)證簡化樹冠阻力模型的可行性，選取4種不同的樹枝進(jìn)行風(fēng)洞實(shí)驗(yàn)以驗(yàn)證數(shù)值模擬結(jié)果。風(fēng)洞實(shí)驗(yàn)裝置如圖2所示，圖中和為壓差測試點(diǎn)。圖3所示為桂花、竹子、廣玉蘭和八角金盤4種典型實(shí)驗(yàn)用樹枝結(jié)構(gòu)，這4種樹種植面廣且形態(tài)差異較大，故將它們作為實(shí)驗(yàn)樹種具有一定代表性。

分別將上述4種樹枝置于測試段(圓形風(fēng)筒)，每種樹枝設(shè)定5種不同葉面積指數(shù)，共計(jì)進(jìn)行20組實(shí)驗(yàn)。葉面積指數(shù)通過控制實(shí)驗(yàn)樹枝的數(shù)量以及葉片數(shù)來實(shí)現(xiàn)。葉面積指數(shù)取值范圍為0~20。根據(jù)葉面積指數(shù)定義可知，測試段7內(nèi)樹枝所有葉片總面積與測試段7底面積之比即為葉面積指數(shù)。實(shí)驗(yàn)時(shí)，通過無級(jí)變頻裝置2(0~50 Hz內(nèi)無級(jí)調(diào)節(jié))和蜂窩柵極6均勻向測試段輸送速度1.0~15.4 m/s的連續(xù)風(fēng)速。并在測試段7的和兩點(diǎn)采用數(shù)字壓差計(jì)測量氣流通過樹枝前后的壓力損失，和兩點(diǎn)分別在測試段圓周上均勻取6個(gè)測試點(diǎn)以測出和兩點(diǎn)平均靜壓差。測試段7尾部通向大氣，每變換一次風(fēng)速等待片刻，待風(fēng)速均勻穩(wěn)定后再進(jìn)行測試。風(fēng)速測量可以通過孔板流量計(jì)間接測得，也可以直接在測試段7使用熱線風(fēng)速儀測量。為獲得更高的測量精度，本實(shí)驗(yàn)采用后者。

1—風(fēng)機(jī)；2—變頻裝置；3—軟接頭；4—孔板流量計(jì)；5—擴(kuò)壓器；6—蜂窩柵極；7—測試段。

圖2 風(fēng)動(dòng)實(shí)驗(yàn)裝置

Fig. 2 Apparatus of wind tunnel for airflow experiment

(a) 桂花；(b) 竹子；(c) 廣玉蘭；(d) 八角金盤

圖3 4種實(shí)驗(yàn)樹枝

Fig. 3 Four kinds of experimental tree branches

3 結(jié)果與討論

為定量描述樹冠的形態(tài)結(jié)構(gòu)，先給出描述樹冠形態(tài)結(jié)構(gòu)的3個(gè)基本參量(樹冠稠密度VF、樹冠分形維數(shù)D、樹冠葉面積指數(shù)AI)的定義及計(jì)算方法。在本文所討論的二維樹冠簡化模型條件下(見圖1)，二維葉面積指數(shù)定義為：所有葉團(tuán)簇周長比二維樹冠底寬，即AI×π×I/0(其中，I為樹葉團(tuán)簇的當(dāng)量直徑，0為二維樹冠底寬，為樹葉團(tuán)簇的數(shù)量)。稠密度VF定義為：所有葉團(tuán)簇的面積比二維樹冠面積，計(jì)算式為VF=×π×I2/(4)(其中為二維樹冠面積，即樹冠輪廓所圍成的面積)。保持樹冠輪廓大小和形狀不變，通過改變?nèi)~團(tuán)簇的當(dāng)量直徑I或樹葉團(tuán)簇的數(shù)量，構(gòu)成不同的樹冠模型，依據(jù)上述給出的樹冠形態(tài)描述參量的計(jì)算式，分別計(jì)算出相應(yīng)模型的樹冠分形維數(shù)、葉面積指數(shù)及稠密度。本文中的樹冠分形維數(shù)f采用盒維數(shù)方法計(jì)算[27]，最后應(yīng)用計(jì)算流體力學(xué)軟件Fluent對(duì)樹冠模型進(jìn)行數(shù)值計(jì)算，討論樹冠阻力系數(shù)與樹冠形態(tài)結(jié)構(gòu)的關(guān)系，并將樹冠流場模擬結(jié)果與文獻(xiàn)結(jié)果進(jìn)行比較，以驗(yàn)證這種樹冠簡化方法的可行性。最后在數(shù)值模擬和實(shí)驗(yàn)結(jié)果基礎(chǔ)上擬合給出樹冠阻力系數(shù)與這3個(gè)基本參量的關(guān)系式。

根據(jù)上述模型與方法以及大量學(xué)者的研究可知：樹冠阻力系數(shù)Cn與樹冠稠密度VF、樹冠分形維數(shù)f以及葉面積指數(shù)AI有關(guān)，即n=(VF，f，AI)。如MOLINA-AIZ等[3]通過風(fēng)洞實(shí)驗(yàn)探究了樹冠阻力系數(shù)與葉面積密度的關(guān)系，GRANT等[28]進(jìn)行現(xiàn)場測試探究了雪松樹冠阻力系數(shù)與雷諾數(shù)以及孔隙率的關(guān)系，而BEREZOVSKAVA等[27]通過實(shí)驗(yàn)、理論以及計(jì)算模擬方法分析了樹冠與分形的關(guān)系。下面分析VF，f與AI之間的關(guān)系，然后選擇與Cn最相關(guān)的無量綱變量。通過分析圖4發(fā)現(xiàn)：簡化模型葉面積指數(shù)AI與樹冠稠密度VF成正比，分形維數(shù)f與樹冠稠密度VF以及葉面積指數(shù)AI都成對(duì)數(shù)關(guān)系，本文采用葉面積指數(shù)AI分析樹冠阻力系數(shù)n，即n=(AI)。

3.1 數(shù)值模型可靠性驗(yàn)證

為驗(yàn)證數(shù)值模型的可靠性，將本文的數(shù)值計(jì)算結(jié)果與其他文獻(xiàn)結(jié)果進(jìn)行比較。圖5(a)給出了葉面積指數(shù)AI=8.3(與文獻(xiàn)中取值一致)情形下，流域中高度為/=1處方向的速度分布。作為對(duì)比，圖5(a)中還給出GREEN[15]和KENJERE?等[29]的實(shí)驗(yàn)結(jié)果。由圖5(a)可知：空氣到達(dá)樹冠迎風(fēng)邊界前，沿方向的流速幾乎保持與入口處流速相同，在/=0處空氣流速達(dá)到最大值，隨后流速迅速降低，在/=4.5處流速出現(xiàn)低谷點(diǎn)；此后流速開始回升，并最終保持在某一速度。與其他研究者的結(jié)果對(duì)比可知：本文給出的樹冠迎風(fēng)區(qū)域的風(fēng)速略高于文獻(xiàn)[15]和[29]中的實(shí)驗(yàn)風(fēng)速，但在樹冠的下游區(qū)域本文數(shù)值結(jié)果與實(shí)驗(yàn)結(jié)果吻合較好。圖5(b)給出了/=1處(樹冠下游)，空氣沿方向的平均速度分布。為方便比較，圖5(b)中橫坐標(biāo)的速度參量采用其他文獻(xiàn)給出的量綱一速度形式，即將實(shí)際風(fēng)速比摩擦速度*(=1.03 m/s)。摩擦速度定義如下：

(a)VF與AI；(b)VF與f；(c)AI與f

圖4 稠密度VF、分形維數(shù)f與葉面積指數(shù)AI相互關(guān)系

Fig. 4 Relationship among solid volume fractionVF, fractal dimensionfand leaf area indexAI

(a)方向速度分布驗(yàn)證(/=1)；(b)方向速度分布驗(yàn)證(/=1)

圖5方向速度分布驗(yàn)證(/=1)和方向速度分布驗(yàn)證(/=1)

Fig. 5 Normalized horizontal profiles of streamwise velocity component (/=1) and normalized vertical profiles of streamwise velocity component (/=1)

式中：為馮卡門常數(shù)(取0.4)；top為樹冠頂部速度；top和0分別為樹冠高度和參考高度?？梢钥闯觯罕疚牡臄?shù)值結(jié)果與文獻(xiàn)[16]中的數(shù)值結(jié)果吻合較好，而文獻(xiàn)[4]和[30]中的計(jì)算結(jié)果均高于本文結(jié)果。引起這些誤差的原因可能有2個(gè)：一是?湍流模型不能詳盡模擬樹冠周圍劇烈變化的湍流分布；二是二維簡化樹冠與三維真實(shí)樹冠在模擬樹冠內(nèi)外流場分布的流通截面積上存在偏差。

上述對(duì)比分析結(jié)果表明：將三維樹冠簡化為二維樹冠后，并采用標(biāo)準(zhǔn)?模型進(jìn)行樹冠流場模擬雖然存在一定的誤差，但仍可獲得較為滿意的結(jié)果。

3.2 樹冠流場特性

由于樹冠枝葉的阻礙作用，空氣流速在樹冠內(nèi)部將降低，伴有明顯的速度梯度和壓力梯度，并產(chǎn)生劇烈的湍流脈動(dòng)。選取模擬結(jié)果中葉團(tuán)簇的當(dāng)量直徑I=4.84 cm，速度入口=1 m/s時(shí)葉面積指數(shù)分別為5.3，8.3，12.7和18.5的4種情形，并與文獻(xiàn)[5]和[24]中的模擬結(jié)果進(jìn)行比較，見圖6。IRVINE等[5]的壓強(qiáng)分布圖表明：來流遇到樹冠時(shí)靜壓由0 Pa迅速升高為3.8 Pa，離開樹冠時(shí)壓強(qiáng)由?1.9 Pa迅速回升到0 Pa；GASH[24]的靜壓分布圖顯示風(fēng)穿過森林時(shí)，森林尾部的靜壓分布為?4.3~0 Pa。而本文模擬結(jié)果與文獻(xiàn)[5]和[24]中的模擬結(jié)果相比：空氣流過樹冠前后壓力梯度分布相似，即氣流在遇到樹冠時(shí)由于樹冠的阻力作用風(fēng)速劇烈衰減，靜壓升高，并在樹冠內(nèi)部壓靜壓發(fā)生轉(zhuǎn)折迅速由正值降為負(fù)值，離開樹冠時(shí)靜壓迅速回升，但是數(shù)值上有一定的差距，靜壓分布為0~0.6 Pa及?0.6~0 Pa。與文獻(xiàn)[5]和[24]中模擬的二維樹冠壓力分布相比，本文模擬值在絕對(duì)值上偏小的原因有： 1) 二維樹冠簡化方式不同，模擬出來的壓力分布也不同，在這點(diǎn)上文獻(xiàn)[5]和[24]中模擬的尾流壓力分布也存在很大差異；2) 文獻(xiàn)[5]和[24]中模擬的是二維簡化森林邊界的流場分布，來流穿過森林前后邊界的壓力分布和來流穿過單個(gè)二維簡化樹冠肯定是有差別的。此外，由圖6可知：在入口風(fēng)速相同的情況下，隨著樹冠葉面積密度增加，靜壓變化的梯度越加明顯。

(a)AI=5.3；(b)AI=8.3；(c)AI=12.7；(d)AI=18.5；(e) 文獻(xiàn)[5]的模擬結(jié)果；(f) 文獻(xiàn)[24]的模擬結(jié)果

圖6 葉團(tuán)簇的當(dāng)量直徑I=4.84 cm，入口速度=1 m/s時(shí)的壓力云圖

Fig. 6 Pressure contour with diameter of leaf clustersI=4.84 cm and inlet velocity=1 m/s

3.3 樹冠阻力系數(shù)n模型

根據(jù)SANZ[31]提出樹冠動(dòng)量源項(xiàng)u：

式中：d為阻力系數(shù)；為空氣密度；f為葉面積指數(shù)；為空氣平均流速。

樹冠阻力的表達(dá)式為

由此提出樹冠阻力系數(shù)n：

由于落葉樹秋季落葉后形成的樹干也會(huì)產(chǎn)生流動(dòng)阻力，故這里包含常數(shù)項(xiàng)，以修正式(8)。將所有的模擬結(jié)果擬合得出：

代入式(9)得到：

由式(11)可知：葉面積指數(shù)AI的二次項(xiàng)系數(shù)相比常數(shù)項(xiàng)和一次項(xiàng)可以忽略，樹冠阻力系數(shù)Cn為葉面積指數(shù)AI的一次函數(shù)，即

本文采用4種植物(桂花、竹子、廣玉蘭和八角金盤)進(jìn)行實(shí)驗(yàn)，每種植物分別取5種葉面積指數(shù)進(jìn)行實(shí)驗(yàn)，則可得到20條不同樹種不同葉面積指數(shù)的壓力損失?速度曲線。分別將這20條曲線所對(duì)應(yīng)的葉面積指數(shù)AI以及速度代入式(13)得出20條模擬曲線。但模擬公式得出的壓差偏小，故須乘以修正系數(shù)，即

這樣得出的模擬結(jié)果與實(shí)驗(yàn)結(jié)果進(jìn)行對(duì)比，結(jié)果如圖7所示。由圖7可知：當(dāng)模擬公式乘以修正系數(shù)后得出的結(jié)果與實(shí)驗(yàn)結(jié)果十分吻合。其中經(jīng)過修正的擬合公式代入各樹種的參數(shù)后和風(fēng)洞試驗(yàn)所得結(jié)果誤差見表2。則樹冠阻力系數(shù)n可修正為′n：

(a) 桂花；(b) 竹子；(c) 廣玉蘭；(d) 八角金盤

圖7 實(shí)驗(yàn)數(shù)據(jù)與經(jīng)過修正的擬合結(jié)果對(duì)比

Fig. 7 Comparison of modified simulation results and measured results

表2 經(jīng)過修正的擬合公式所得結(jié)果與實(shí)驗(yàn)結(jié)果誤差評(píng)估

4 結(jié)論

1) 為深入分析樹冠內(nèi)部結(jié)構(gòu)參數(shù)對(duì)樹冠繞流流動(dòng)阻力的影響，將三維復(fù)雜樹冠簡化為二維簡化樹冠模型并提出相關(guān)樹冠結(jié)構(gòu)參數(shù)葉面積指數(shù)AI、稠密度VF以及分形維數(shù)f。采用數(shù)值模擬和實(shí)驗(yàn)測試相結(jié)合的方法研究幾種典型形態(tài)特征樹冠的內(nèi)部流場特征和流動(dòng)阻力。用樹冠葉面積指數(shù)(AI)、稠密度(VF)和分形維數(shù)(f)表征樹冠復(fù)雜的形態(tài)結(jié)構(gòu)，討論了樹冠形態(tài)結(jié)構(gòu)參數(shù)對(duì)樹冠流動(dòng)阻力的影響規(guī)律。

2) 通過與文獻(xiàn)值對(duì)比，可以驗(yàn)證簡化樹冠水平方向(方向)速度分布與文獻(xiàn)中的實(shí)驗(yàn)速度及模擬速度較吻合；簡化樹冠高度方向(方向)的平均速度分布曲線與其他模擬結(jié)果及文獻(xiàn)中的實(shí)驗(yàn)結(jié)果相比一致 性高。

3) 在簡化樹冠的基礎(chǔ)上提出樹冠阻力模型Δ=n×2/2，并假設(shè)樹冠阻力系數(shù)n是模型樹冠結(jié)構(gòu)參數(shù)稠密度VF、分形維數(shù)f及葉面積指數(shù)AI的函數(shù)，即n=f(VF，f，AI)，通過分析發(fā)現(xiàn)簡化模型葉面積指數(shù)AI與樹冠稠密度VF成正比，分形維數(shù)f與樹冠稠密度VF以及葉面積指數(shù)AI都成對(duì)數(shù) 關(guān)系。

4) 采用葉面積指數(shù)AI分析樹冠阻力系數(shù)n，進(jìn)而引出樹冠阻力系數(shù)模型n=+×AI+×AI2。通過將模擬結(jié)果擬合分析得出關(guān)于樹冠阻力系數(shù)n的模型公式：′n=×C=×(0.05+0.012AI)。

[1] CAO Jinxin, TAMURA Y, YOSHIDA A. Wind tunnel study on aerodynamic characteristics of shrubby specimens of three tree species[J]. Urban Forestry & Urban Greening, 2012, 11: 465?476.

[2] SANZ R J, van BEECK J, DEZSO W G. Wind tunnel simulation of the wind conditions inside bidimensional forest clear-cuts. Application to wind turbine siting[J]. Journal of Wind Engineering and Industrial Aerodynamics, 2007, 95: 609?634.

[3] MOLINA-AIZ F D, VALERA D L, ALVAREZ A J, et al. A wind tunnel study of airflow through horticultural crops: determination of the drag coefficient[J]. Biosystems Engineering, 2006, 93(4): 447?457.

[4] DUMAN T, KATUL G G, SIQUEIRA M B, et al. A velocity–dissipation lagrangian stochastic model for turbulent dispersion in atmospheric boundary-layer and canopy flows[J]. Boundary Layer Meteorology, 2014, 152: 1?18.

[5] IRVINE M R, GARDINER B A, HILL M K. The evolution of turbulence across a forest edge[J]. Boundary Layer Meteorology, 1997, 84(3): 467?496.

[6] NOVAK M D, WARLAND J S, ORCHANSKY A L, et al. Wind tunnel and field measurements of turbulent flow in forests. Part I: uniformly thinned stands[J]. Boundary Layer Meteorology, 2000, 95(3): 457?495.

[7] 關(guān)德新, 朱廷耀. 樹冠結(jié)構(gòu)參數(shù)及附近風(fēng)場特征的風(fēng)洞模擬研究[J]. 運(yùn)用生態(tài)學(xué)報(bào), 2000, 11(2): 202?204. GUAN Dexin, ZHU Tingyao. Wind tunnel experiment on canopy structural parameters of isolated tree and wind velocity field characters nearby[J]. Chinese Journal of Applied Ecology, 2000, 11(2): 202?204.

[8] ZENG P, TAKAHASHI H. A first-order closure model for the wind flow within and above vegetation canopies[J].ultural andogy, 2000, 103: 301?313.

[9] SIQUEIRA M, KATUL G, TANNY J. The effect of the screen on the mass, momentum, and energy exchange rates of a uniform crop situated in an extensive screen house[J]. Boundary Layer Meteorology, 2012, 142(3): 339?363.

[10] SOGACHEV A, PANFEROV O. Modification of two-equation models to account for plant drag[J]. Boundary Layer Meteorology, 2006, 121(2): 229?266.

[11] SOGACHEV A. A note on two-equation closure modelling of canopy flow[J]. Boundary Layer Meteorology, 2009, 130(3): 423?435.

[12] MOCHIDA A, TABATA Y, IWATA T, et al. Examining tree canopy models for CFD prediction of wind environment at pedestrian level[J]. Journal of Wind Engineering and Industrial Aerodynamics, 2008, 96: 1667?1677.

[13] HIRAOKA H, OHASHI M. A (?) turbulence closure model for plant canopy flows[J]. Journal of Wind Engineering and Industrial Aerodynamics, 2008, 96: 2139?2149.

[14] 李亮, 李曉峰, 林波榮, 等. 用帶?源項(xiàng)兩方程湍流模型模擬樹冠流[J]. 清華大學(xué)學(xué)報(bào)(自然科學(xué)版), 2006, 46(6): 753?756. LI Liang, LI Xiaofeng, LIN Borong, et al. Simulation of canopy flows using?two-equation turbulence model with source/sink terms[J]. Journal of Tsinghua University (Science and Technology), 2006, 46(6): 753?756.

[15] GREEN S R. Modeling turbulent air flow in a stand of widely spaced trees[J]. Phoenics, 1992, 23: 294?312.

[16] HUANG J, CASSIANI M, ALBERTSON J D. Coherent turbulent structures across a vegetation discontinuity[J]. Boundary Layer Meteorology, 2011, 140: 1?22.

[17] DUPONT S, BRUNET Y. Edge flow and canopy structure: a large-eddy simulation study[J]. Boundary Layer Meteorology, 2008, 126: 51?71.

[18] DUPONT S, BRUNET Y. Influence of foliar density profile on canopy flow: a large-eddy simulation study[J]. Agricultural and Forest Meteorology, 2008, 148: 976?990.

[19] PATTON E G, KATUL G G. Turbulent pressure and velocity perturbations induced by gentle hills covered with sparse and dense canopies[J]. Boundary Layer Meteorology, 2009, 133: 189?217.

[20] YANG B, RAUPACH M R, SHAW R H, et al. Large-eddy simulation of turbulent flow across a forest edge. Part I: flow statistics[J]. Boundary Layer Meteorology, 2006, 120(3): 377?412.

[21] SCHLEGEL F, STILLER J, BIENERT A, et al. Large-eddy simulation of inhomogeneous canopy flows using high resolution terrestrial laser scanning data[J]. Boundary Layer Meteorology, 2012, 142(2): 223?243.

[22] BELCHER S E, HARMAN I N, FINNIGAN J J. The wind in the willows: flows in forest canopies in complex terrain[J]. Annual Review of Fluid Mechanics, 2012, 44(1): 479?504.

[23] SVENSSON U, HAGGKVIST K. A two-equation turbulence model foor canopy flows[J]. Journal of Wind Engineering and Industrial Aerodynamics, 1990, 35: 201?211.

[24] GASH J H C. Observations of turbulence downwind of a forest-heath interface[J]. Boundary Layer Meteorology, 1986, 36(3): 227?237.

[25] LAUNDER B E, SPALDING D B. The numerical computation of turbulent flows[J]. Computer Methods in Applied Mechanics and Engineering, 1974, 3: 269?289.

[26] VERSTEEG H K, MALALASEKERA W. An introduction to computational fluid dynamics[M]. Englewood Cliffs, NJ: Prentice-Hall, 1995: 82.

[27] BEREZOVSKAVA F S, KAREV G P, KHLEBOPROS R G, et al. A fractal approach to computer analytical modelling of tree crowns[J]. Trees, 1997, 11: 323?327.

[28] GRANT P F, NICKLING W G. Direct field measurement of wind drag on vegetation for application to windbreak design and modelling[J]. Land Degradation & Development, 1998, 9: 57?66.

[29] KENJERE? S, KUILE B T. Modelling and simulations of turbulent flows in urban areas with vegetation[J]. Journal of Wind Engineering and Industrial Aerodynamics, 2013, 123: 43?45.

[30] REYNOLDS A M. Development and validation of a lagrangian probability density function model of horizontally-homogeneous Turbulence within and above plant canopies[J]. Boundary Layer Meteorology, 2012, 142: 193?205.

[31] SANZ C. A note on?modelling of vegetation canopy air-flows[J]. Boundary Layer Meteorology, 2003, 108: 191?197.

(編輯 楊幼平)

Numerical simulation and experimentalstudy of canopyflowresistance characteristics

YANG Hui1, 2, FU Haiming1

(1. School of Environmental Science and Engineering, Donghua University, Shanghai 200051, China;2. Department of Building Environment and Energy Engineering,Guilin University of Aerospace Technology, Guilin 541004, China)

Computational fluid dynamics (CFD) and field experiments were used to investigate the flow characteristics and flow resistance through vegetation canopies with several typicalmorphologicalcharacteristics. The leaf area index (AI), solid volume fraction (VF) and thefractal dimension(f) were used to analyze the complicated morphology of vegetation canopy. The influences of leaf area index, solid volume fraction and fractal dimension on flow resistance were studied. The results show that leaf area index approximately has linear relationship with solid volume fraction, and thefractal dimension of canopy has logarithmic relationship with leaf area index and solid volume fraction, respectively. A model expression between the drag coefficient of the canopy and leaf area index was presented by analyzing simulation results. The model expression was validated using experimental results in a wind tunnel with four kinds of tree branches. The results show that the simulation results from 2D simplified model are slightly larger than experimental results, and can qualitatively predict the trend of the relationship between the drag coefficient of the canopy and leaf area index. However, the results obtained by using the proposed 2D simplified canopy model show good agreement with the experimental results after modifying the simulation results with a correction factor.

canopy; flow resistance; drag coefficient of canopy; leaf area index

10.11817/j.issn.1672-7207.2016.12.043

X513; S731.2

A

1672?7207(2016)12?4292?09

2015?12?28；

2016?03?09

國家自然科學(xué)基金資助項(xiàng)目(41371445)(Project(41371445) supported by the National Natural Science Foundation of China)

付海明，博士，高級(jí)工程師，從事空氣過濾凈化與空氣品質(zhì)研究；E-mail：fhm@dhu.edu.cn