零方程湍流模型用于城市微氣候的模擬研究

劉吉營，郭敏，斯雷布里奇·伊蓮娜，達維多維奇·丹科

（1.山東建筑大學熱能學院，山東濟南250101；2.山東建筑大學山東省建筑節能技術重點實驗室，山東濟南250101；3.山東省可再生能源建筑應用技術重點實驗室，山東濟南250101；4.美國馬里蘭大學機械工程系，美國馬里蘭州20742；5.胡貝爾工程林木公司，美國商務城30530）

零方程湍流模型用于城市微氣候的模擬研究

劉吉營1，2，3，*，郭敏1，2，3，斯雷布里奇·伊蓮娜4，達維多維奇·丹科5

（1.山東建筑大學熱能學院，山東濟南250101；2.山東建筑大學山東省建筑節能技術重點實驗室，山東濟南250101；3.山東省可再生能源建筑應用技術重點實驗室，山東濟南250101；4.美國馬里蘭大學機械工程系，美國馬里蘭州20742；5.胡貝爾工程林木公司，美國商務城30530）

室外微氣候環境對城鎮人居環境和建筑能耗有重要影響，數值模擬技術作為快速的預測方法能有效評價社區規劃的優劣。文章選用風洞實驗數據對湍流渦粘度表達式進行優化，提出新的零方程湍流模型（Zero-EQ），并與兩方程湍流模型（MMK）、大渦模擬（LES）和風洞實驗數據進行對比分析；通過對壽安新城區零方程湍流模型的應用分析，探討該模型應用于實際工程項目的可行性。結果表明:在同等計算機配置下，Zero-EQ模型比MMK模型快60%，是LES模型計算所需時間的1／15，能有效地預測不規則建筑群的流場分布；通過分析不同地形區域的風熱環境特點，Zero-EQ模型可快速預測城鎮微氣候環境，為城區規劃設計提供合理建議。

零方程湍流模型；城市微氣候環境；計算流體動力學；城市規劃設計

Key words:zero-equation turbulence model；microclimate urban environment；computational fluid dynamics；urban planning and design

0 引言

在城鎮化快速發展的過程中，室外人居的風、光、熱、污染物等微氣候生態環境日益受到廣泛關注，而數值模擬技術作為最有效的預測方法之一，已成為評價低碳生態社區規劃優劣的重要方法［1-3］。目前，應用于城市微氣候環境數值模擬最廣泛的兩種湍流模型分別是基于雷諾平均的 RANS（Reynolds-Average Navier-Stokes Equations）模型和大渦模擬LES（Large Eddy Simulation）。時至今日，由于計算機資源的限制，LES仍限于簡單的建筑模型或者相對較小的區域。通過對比分析常用湍流模型預測城市風環境的優劣，在較小CPU能耗情況下，RANS模型在非近壁區和非穩態流動時模擬結果不理想，而LES需要較大的CPU消耗才能取得相對較好的預測結果［4］。而在實際情況下，對于城鎮區域大小的空間尺度，選擇RANS模型仍需要較大的計算機資源才能獲得可以接受的預測結果。因此，如何能夠快速且相對準確地預測城市微氣候環境狀況，提供給工程師和設計師準確的信息，滿足工程需要，將成為本課題的重點研究方向。近年來，有一種已經被證明能夠快速有效地模擬室內微氣候環境的零方程湍流模型（Zero-EQ），正逐步應用于室外微氣候環境的預測［5-8］。

目前，已有的Zero-EQ模型大多利用普朗特混合長度理論，而混合長度理論主要應用于自由剪切流，適用于簡單的邊界層流動。因此，Van Driest提出了壁面衰減函數來解決管道流的流動問題［9］，Cebeci等提出了一個綜合考慮高速流動、低雷諾流動和管道流等兩層的Zero-EQ模型［10］，而Baldwin等則提出了類似的兩層Zero-EQ模型來解決管道流［11］。Chen等通過設置特定的湍流動能，提出了新的基于室內環境研究的渦粘度表達式［5］。以上Zero-EQ僅適用于特定的邊界層流動或小尺度下的室內環境的研究，對于城鎮尺度下的大氣邊界層，顯然不適用。現有應用于室外微氣候環境預測的Zero-EQ是基于單棟建筑的風洞試驗推導出來的，僅局限于室外單棟建筑或規則建筑群的微氣候環境［8］。對于不同高低、間距和位置等不規則分布的建筑群或實際城鎮社區，因其復雜的建筑擾流與多變的特征長度尺度，目前鮮有涉及［12］。因此，文章通過分析風洞實驗數據，研究不規則建筑群的流場變化，推導出優化的湍流渦粘度表達式，提出新的零方程湍流模型。通過與兩方程湍流模型MMK（Murakami-Mochida-Kondo）、LES和風洞實驗對比，驗證其具有相對準確和快速的優勢，并應用于某城鎮低碳生態建設規劃設計，具有非常廣泛的應用前景［13］。

1 零方程模型優化

基于非穩態不可壓連續性和動量的雷諾平均Navier-Stokes方程（RANS）可用式（1）或式（2）表示為

式中:ui為分量速度，m／s；xi為分量方向，m；t為時間，s；ρ為流體密度，kg／m3；μt為湍流渦粘性系數，Pa·s；β為熱膨脹系數，1／K；為雷諾應力，Pa；T0為操作溫度，℃；T為流體溫度，℃；gi為重力加速度，m／s2。

在式（1）（2）中，由于雷諾應力是一個未知參量，為了能解上述湍流方程，雷諾應力選用Boussinesq假設，把雷諾應力與平均應變率和湍流渦粘度（μt）相關聯（即常說的線性渦粘度模型），得到式（3）:

式中:δij是克羅內克符號；k是湍流動能，m2／s2；ρδijk是一種用張量代數來表達的線性本構關系，用于求解雷諾應力。

為了能夠求解雷諾應力，對兩方程湍流模型而言，湍流渦粘度（μt）與湍流動能（k）和湍流耗散率（ε）或者湍流特定耗散率（ω）相關。因此，為了同時解k和ε，需要兩個封閉的k和ε方程來計算，而對于兩者不同的處理策略，形成了目前眾多不同的兩方程湍流模型。如果只選用k方程來處理雷諾應力，即為一方程湍流模型；如果不選任何方程，只是采用代數算法公式來表述雷諾應力的話，即稱算法模型，也就是Zero-EQ模型。

Zero-EQ模型在利用Boussinesq假設來求解雷諾應力的過程中，忽略k微分方程，而采用簡單算法方程，則式（3）簡化為式（4）

為求解湍流渦粘度（μt），依據Prandtl的混合長度理論，將式（4）簡化為單向的二維邊界層，其平均的雷諾應力與平均速度梯度有關。在本實驗中，僅考慮雷諾在垂直方向上的順流速度梯度，式（4）可進一步簡化為式（5）

式中:U為垂直方向上的順流速度分量，m／s。通過處理建筑周圍流場各測點數據，得到相應的雷諾應力和順流速度分量梯度，計算出各測點的湍流渦粘度分布。通過分析渦粘度數據，推導出優化的零方程湍流渦粘度表達式，如式（6）至式（8）所示

式中:Returb是湍流雷諾數；Rebulk是主導雷諾數；L是積分湍流尺寸長度，m；H為區域內平均建筑高度，m；V是建筑周圍流場局部計算點的速度大小，m／s；kin，H是在建筑H高度處入口湍流動能，m2／s2；τin，H是在建筑高度H處入口湍流時間尺度，s；Uin，H是建筑高度H處入口速度大小，m／s；v是動力粘度，m2／s；a和b為渦粘度表達式系數。

2 湍流模型對比分析

2.1 湍流模型介紹

當標準k-ε湍流方程預測流體經過鈍體建筑時，存在一些誤差，比如過高預測碰撞區域的湍流動能k等，因此Kato等通過修改湍流動能產出項Pk，提出修正的兩方程湍流模型（LK）［14］；但是LK模型在研究湍流雷諾應力和Pk時仍存在數值的不穩定性現象。因此，Tsuchiya等通過在渦粘度表達式中添加修正參數，而非湍流動能產出項，得到了基于LK模型修正的MMK模型［13］。考慮到LES計算需要大量的計算機資源，作為修正的k-ε湍流模型:LK和MMK模型，因其相對快速而又相對準確地預測建筑環境，目前已被大量地應用于城市風熱環境的研究中［15］。在LES計算中，選擇Smagorinsky常數Cs＝0.12用于求解亞格子渦粘度。文章選用MMK模型、LES模型，與Zero-EQ模型和風洞實驗數據對比，分析Zero-EQ模型的優劣。

2.2 風洞實驗介紹

文章選用賓夕法尼亞州立大學航空航天工程系的風洞實驗室開展實驗［12，16-17］（如圖1所示），風洞試驗臺長度約6 m、寬度為0.6 m、高度為0.9 m。在Qian研究單棟建筑的基礎上［8］，選擇四棟毗鄰的學生公寓作為研究對象，其中每棟建筑長28 m、寬20m和高30 m，建筑模型縮尺比例為1:250，雷諾數為80000。

縮尺模型周圍的水平和豎直測點分布圖如圖2所示。在選擇水平測點時，盡量選擇建筑中心側的位置。根據縮尺模型分布和周圍流場的特點，布置了水平位置上的八橫行（Y1至Y8）和八縱列（X1至X8）。研究選取對模型下游流場影響比較大的四橫行Y3、Y4、Y5和Y6作為研究對象，除因位置原因無法測試外，總計26個水平測點。其中測點X1Y1、X1Y4、X1Y8主要用于入口邊界的速度和湍流強度邊界條件。

2.3 CFD模擬設定

在數值計算時，減少計算區域上游的距離滿足實際風洞邊界測點與第一個建筑的距離要求，即計算區域距上游第一個建筑1 H長度，距下游最后一個建筑12 H長度，而側面和垂直區域受實際風洞實驗條件的限制，分別為2.5 H高度和3.6 H高度，計算模型如圖3所示。MMK模型和Zero-EQ模型結構化網格數為140×110×100。采用網格收斂指數GCI對網格相關性分析，即:GCI［fine］＝FS｜ε｜／（rp-1）。計算上游第一個建筑迎風面的表面摩擦系數Cf，分別選取粗糙網格（120×90×85，Δ1），中等密度網格（144×110×100，Δ2）和加密網格（170×132×120，Δ3）。結果表明:GCI［3，2］＝3.18%（加密網格對比中等密度網格）和GCI［2，1］＝7.54%，GCI［3，2］相比GCI［2，1］有較小的網格相關性；因此本模型選取中等網格密度。入口邊界平均流速和湍流強度拋物線如圖4所示。出口采用出口（OUTLET）類型，區域側面采用壁面（WALL）邊界條件，計算區域頂部采用對稱邊界（Symmetry）條件。

圖1 風洞試驗圖

圖2 不同位置測點分布圖

采用LES模型時，絕大部分近壁區滿足y+≈1的要求，網格數為300×200×120，入口邊界采用高斯概率密度分布公式來生成瞬時脈動速度，如式（9）所示

式中:σ為在入口邊界不同高度處的湍流強度，%。時間步長設定為0.0005 s，滿足絕大部分計算區域Courant-Friedrichs-Levy（CFL）值低于1的要求。

圖3 計算模型圖

圖4 入口邊界沿建筑高度方向變化的速度和湍流強度分布圖（a）速度；（b）湍流強度

2.4 計算結果對比分析

在相同計算工況和同等計算機配置條件下，三種湍流模型所消耗的計算時間比為:Zero-EQ: MMK:LES＝1:1.6:15，很顯然，Zero-EQ能夠比MMK湍流模型快60%，是LES模型計算所需時間的1／15。計算三個湍流模型預測結果與風洞實驗數據的平均相對誤差，分別得到12.9%、7.9%和11.8%（MMK，LES和Zero-EQ），可以看出Zero-EQ模型具有一定的優勢。因此，Zero-EQ模型能夠比較快速地模擬建筑群周圍的流場分布，同時Zero-EQ模型能相對準確地預測流場速度大小分布。

文章僅考慮了Y4橫行中的四個測點X2Y4，X4Y4，X5Y4和X8Y4的縱向分速度分布，更多橫行對比分析可參考文獻［17］。因Y4行位于多棟建筑的中心位置，該行測點的對比分析能夠有效地評估Zero-EQ模型的優劣。圖5為各測點的MMK、LES、Zero-EQ的模擬結果與風洞實驗數據對比，在縱向分速度上，三個模型皆能準確地預測測點X2Y4的速度，但Zero-EQ能更準確地預測建筑頂部區域的漩渦。在測點X4Y4，MMK和LES雖能預測流場分布的趨勢，但其速度大小較風洞實驗數值偏小，而Zero-EQ雖與實驗數據誤差較小，卻不能有效捕捉流場分布；在測點X5Y4，三個模型的預測結果與風洞實驗數據吻合較好；在測點X8Y4，MMK和LES模型的預測速度與實驗結果相比偏小，而Zero-EQ雖不能完全地預測出建筑背風面區域流場的渦旋分布，但能有效地預測此區域的速度大小。

從風洞實驗數據及MMK和LES模型預測結果對比，可以看出湍流模型能相對準確地預測流場速度大小的分布，但不能完全地預測建筑擾流場的流動趨勢，這是由于（1）現有的湍流模型仍不能準確地描述非近壁區的轉捩流動以及非穩態流動中的漩渦脫落；（2）風洞試驗采用的熱線風速儀不能完全準確地捕捉非穩態流動漩渦區域的復雜流動，這是試驗的不足之處，今后將采用更精確的實驗設備。但總的說來，通過驗證和分析，相比MMK和LES模型，Zero-EQ能夠快速而又相對有效地模擬不規則建筑群的流場分布。

3 零方程湍流模型應用研究

3.1 低碳生態建設專項規劃介紹

文章選擇成都市蒲江縣壽安鎮生態建設項目，來探討零方程湍流模型應用于實際工程項目的可行性。結合壽安鎮當地氣象條件，建立較準確的城市三維模型，分析城區各區域風熱環境的特點，提出壽安新城的城市規劃建議，達到工程應用的目的。

雖然壽安鎮舊城區作為一個小區域，但其仍滿足城市地形區域UTZ（Urban Terrain Zone）的部分分類特點［18］，即:鎮中心商業區（UTZ-1）、壽安中學校區（UTZ-2）、高層樓居民生活區（UTZ-5）和鎮郊居民區（UTZ-6）。其各自詳細的建筑信息見表1和如圖6所示。

圖5 不同測點縱向速度（U）沿建筑高度方向（Z坐標）變化的分布圖

表1 壽安舊城區不同區域的特點

3.2 模擬設定

根據近30多年浦江縣的氣象資料，年平均風速為0.9 m／s，年靜風頻率超過40%，空氣水平擴散條件差，全年主導風向為NNE，年平均氣溫為16.5℃。建立壽安鎮準確的城區建筑模型，分別設定計算區域入口邊界、太陽輻射模型和壁面邊界條件:“WIND”入口條件，考慮太陽輻射對地面的影響，Tin＝25℃，選擇夏季6月23日正午時分，入口速度選擇對數函數分布輪廓線，U10＝0.9 m／s；建筑群表面設置溫度邊界條件T＝30℃，表面粗糙度為0.005 m；計算區域網格為600×500×120，滿足絕大部分壁面無量綱數y+＜300。

3.3 結果分析

根據四個不同地形區域建筑分布的特點，分析區域內風熱環境的分布，可供建筑師和工程師做規劃設計時參考。區域內速度大小和溫度分布圖分別如圖7、8所示，分布圖高度選擇行人的肩膀至頭部位置處，大約為1.5 m。分析如下:

（1）鎮中心商業區（UTZ-1）由于商業區位于城鎮中心，街道普遍較窄（約7 m），受周圍較高建筑的影響大，由圖7可知，區域內的速度值普遍偏低，低于0.20 m／s；而街區內溫度較高，如圖8所示，高于35°C，部分位置溫度甚至超過40°C。很顯然，在街區通風不暢、太陽輻射比較強烈的情況下，街區內的熱環境舒適度比較差，行人普遍感到不舒適；受此影響，城市熱島現象比較嚴重，區域內溫度整體偏高，使得商業區的建筑能耗增大。

圖6 壽安舊城區建筑模型和區域分類圖

（2）壽安中學校區（UTZ-2）在壽安中學校區內，因學生操場位置周圍無高層建筑，操場內風速大小受周圍環境影響較小，與城鎮近郊的風速相近，最高達到0.45 m／s；學生教學樓和教師辦公樓周圍流場風速較鎮中心商業區高，風環境舒適度較好，受其影響，此區域溫度大多維持在28～32°C之間，最高不超過36°C，熱環境舒適度基本可以接受。

圖7 壽安舊城區速度大小分布圖在1.5 m處

圖8 壽安舊城區溫度分布圖在1.5 m處

（3）高層樓居民生活區（UTZ-5）在高層樓居民區，由于其遠離鎮中心，且在主導風向上游區域，風環境比鎮中心商業區稍好，但由于大多數的高層住宅樓選用“口”字型設計，使得幾棟樓中心區域的風速偏低，普遍低于0.1 m／s，風環境舒適度較差；相應地，高層樓居民區內的溫度偏高，維持在35°C以上，最高可達43°C，行人區域的熱環境較差。

（4）鎮郊居民區（UTZ-6）在城鎮近郊，居民房屋大多呈單棟位置分布，房屋之間的間距較大，且房屋高度一般在4～7 m之間，遮蔽效應不明顯，因此居民區周圍的風環境較好（大約0.4 m／s），熱環境舒適度較好，溫度普遍維持在30°C上下。

通過對壽安鎮舊城區不同地形區域風熱環境分析，對新城區城市規劃提供如下建議:

（1）有必要通過合理規劃商業區的布局，采取增大街區寬度，加大綠化、水體建設等措施來減少城市熱島效應的影響；同時，根據當地盛行風向和風速大小，選擇偏南北走向的街道布局，加強城鎮建筑間通風。

（2）受當地風速偏低的影響，嚴格規劃城鎮高層居民區布局設計，杜絕“口”字型分布設計，合理控制建筑群的密度，保證行人區域內風環境舒適性要求；

（3）建議今后在新城規劃設計時，對不同建筑規劃設計方案建立CFD建筑模型，分析其受當地氣象條件影響的微氣候環境，提供給建筑師和工程師足夠準確的信息。

4 結論

通過上述研究可知:

（1）Zero-EQ模型雖不能完全預測建筑擾流的流動趨勢，但能相對準確地預測流場速度分布；在同等計算機配置下，Zero-EQ模型比MMK模型快60%，是LES模型計算所需時間的1／15，能快速且相對有效地預測不規則建筑群的風環境分布。

（2）通過分析不同地形區域的風熱環境特點，Zero-EQ模型能夠快速地預測城鎮微氣候環境，可供建筑師和工程師做城市規劃和生態環境設計時參考，達到實際工程應用的目的。

（3）建議今后城鎮低碳社區在建設過程中，選用CFD模擬技術對新建或改建建筑周圍微氣候環境進行預測，優化和改善城鎮微氣候環境狀況。

文章的零方程湍流模型是基于縮尺風洞試驗得到的，而此風洞試驗未考慮強浮力流動效應等影響，今后將繼續研究零湍流模型中自然對流引起的浮力流動對建筑擾流的影響，同時將對城鎮社區尺度下的微氣候環境進行現場實測，驗證零方程湍流模型的準確性。

［1］ 關吉平，任鵬杰，周成，等.高層建筑行人高度風環境風洞試驗研究［J］.山東建筑大學學報，2010，25（1）:21-25.

［2］ 閆整，房文娟，肖華斌.上海低碳生態城市的規劃實踐進展［J］.山東建筑大學學報，2011，26（3）:266-270.

［3］ 牛盛楠，張欣宜，黃成，等.天津地區居住區采光宇室外風環境模擬研究［J］.山東建筑大學學報，2013，28（1）:12-17.

［4］ Murakami S.，Ooka R.，Mochida A.，et al.CFD analysis of wind climate from human scale to urban scale［J］.Journal ofWind Engineering and Industrial Aerodynamics，1999，81:57-81.

［5］ Chen Q.，Xu W..A zero-equation turbulence model for indoor airflow simulation［J］.Energy Buildings，1998，28（2）:137-144.

［6］ Zhai Z.Q.，Chen Q.Y.，Haves P.，et al.On approaches to couple energy simulation and computational fluid dynamics programs［J］.Build Environment，2002，37（8-9）:857-864.

［7］ Srebric J.，Chen Q.Y.，Glicksman L..Validation of a zeroequation turbulence model for complex indoor airflow simulaiton［J］.ASHRAE Transactions，1999，105（2）:414-427.

［8］ Qian Y..Development of algebraic turbulence model for airflow and contaminant simulations around building［D］.University Park:The Pennsylvania State University，2004.

［9］ Van Driest E.R..On turbulent flow near a wall［J］.Journal of Aero Science，1956，23（1）:1007-1010.

［10］Cebeci T.，Smith A..Analysis of turbulence boundary layers，series in applied mathematics and methods:XV［M］.University Park:Academic Press，1974.

［11］Baldwin B.S.，Lomax H..Thin-layer approximation and algebraic models for separated turbulence flows［J］.AIAA Paper，1978，278-286.

［12］Davidovic D..Improvements in numerical airflow modeling around multiple buidlings［D］.University Park:The Pennsylvania State University，2009.

［13］Tsuchiya M.，Murakami S.，Mochida A.，et al.Development of a new kappa-epsilonmodel for flow and pressure fieldsaround bluff body［J］.Journal of Wind Engineering and Industrial Aerodynamics，1997，67-68（1）:169-182.

［14］Kato M.，Launder B.E..Themodeling of turbulent flow around stationary and vibrating square cylinders［J］.Preprints of 9th Symposium on Turbulent Shear Flow，1993，10（4）:1-6.

［15］Tominaga Y.，Mochida A.，Yoshie R.，et al.AIJguidelines for practical applications of CFD to pedestrian wind environment around buildings［J］.Journal ofWind Engineering and Industrial Aerodynamics，2008，96（10-11）:1749-1761.

［16］Davidovic D.，Liu J.，Heidarinejad M.，etal.Airflow study for a cluster of campus buildings using different turbulence modeling approaches［J］.Invited Paper for the Third Issue of International Journal in Building，Urban，Interior and Landscape Technology，2014.

［17］劉吉營.基于零方程湍流模型的非等溫城市區域熱環境的模擬研究［D］.成都:西南交通大學，2013.

［18］Ellefsen R..Current assessmentof Building Construction Types in Worldwide Example Cities［M］.Dahlgren:Naval Surface Warfare Center Dahlgren Division，1999.

（責任編輯：吳芹）

Numerical simulation of them icro climate urban environment using a zero-equation turbulencemodel

Liu Jiying1，2，3，*，Guo Min1，2，3，Srebric Jelena4，et al.

（1.School of Thermal Engineering，Shandong Jianzhu University，Jinan 250101，China；2.Shandong Key Laboratory of Building Energy Saving Technologies，Shandong Jianzhu University，Jinan 250101，China；3.Shandong Key Laboratory of Renewable Energy Utilization Technologies in Building，Shandong Jianzhu University，Jinan 250101，China；4.Department of Mechanical Engineering，University of Maryland，University Park 20742，USA）

The micro climate urban environment has an important effect on outdoor residential environment and building energy consumption.Numerical simulation as a fast prediction method can effectively evaluate the quality of neighborhood planning.In this paper，a rapid and relatively reliable Zero-equation turbulencemodel（Zero-EQ）is provided through analyzing and improvingmeasurement data of the wind tunnel experiment，and compared with two-equation turbulence model（MMK），Large Eddy Simulation（LES）and a wind tunnel experiment.By investigating the application of Zero-EQmodel in Shouan Town，this research discusses the feasibility of themodel's application to the real project.The result shows that the computational speed of Zero-EQmodel is faster by 60%than MMK model，and only 1／15 of thatusing LESmodel under the same computational resources，indicating the capability of Zero-EQ on relatively accurately predicting the wind flow field in a complex neighborhood.In the end，Zero-EQ can quickly predict the micro climate urban environment and provide reasonable suggestion for the urban planning and design by means of studying the wind and thermal environment over different urban terrains.

TU375.1

A

1673-7644（2014）06-0512-08

2014-10-31

國家科技支撐計劃子課題項目（2012BAJ06B03-01）；山東建筑大學博士科研基金項目（XNBS1408）

劉吉營（1983-），男，講師，博士，主要從事城市熱環境和綠色建筑等方面的研究.E-mail:jxl83＠sdjzu.edu.cn