迷宮流道滴頭內流場和顆粒運動的不同湍流模型數值模擬

唐學林，趙旭紅，李云開，李小芹,2

（1. 中國農業大學水利與土木工程學院，北京 100083；2. 北京市供水管網系統安全與節能工程技術研究中心，北京 100083）

0 引 言

滴灌是目前微灌系統中最有效、最先進的一種節水灌溉方式[1]，滴頭內部流道結構錯綜復雜，消能效果非常好，但也存在滴頭堵塞嚴重、灌水均勻度不高等問題。目前針對滴頭內迷宮流道堵塞這一問題，國內外學者進行了大量的試驗及模擬[2-5]，同時也對迷宮流道的結構作出了優化改進，但堵塞問題依然存在[4]，加上灌區水源多為井水、河水等近地表水，水質較差，使得灌水器堵塞已成為制約滴灌技術應用和推廣的瓶頸[6]。優化滴頭的抗堵塞性能，必須對其流道內部固-液兩相的運動機理有一個較為全面的了解，進而分析滴頭流道內部的流動，實現流道內部流體、顆粒流動的可視化[7]。

由于迷宮流道斷面尺寸微小，一般在0.5～1.2 mm之間，介于微尺度（<1 mm）和常規尺度（>1 mm）之間[8-9]，加上結構復雜，使得迷宮流道滴頭極易堵塞，為了更好地研究灌水器內部的臨界尺度流動問題，在基于宏觀尺度研究方法的CFD湍流模型流動分析的基礎上，有待于依靠試驗測試進行驗證[10]。魏正英等[11-13]用激光多普勒測速儀(LDV)測量了流道中的速度流場，驗證了流態模擬計算的正確性，之后又利用微 PIV對流道內兩相流流場進行試驗，指出造成流道堵塞的主要原因是流道內的漩渦區和低速區。喻黎明等[14-16]利用 PIV觀測迷宮流道的抗堵塞性能進行驗證數值模擬結果的準確性，并提出有效的抗堵塞優化設計方法，但對于致堵顆粒的運動研究多通過數理統計方法分析不同性質顆粒對于流道堵塞的影響。大量專家學者利用不同湍流模型開展迷宮流道內流場水力性能的分析，常用的湍流模型有標準k-ε模型、RNGk-ε模型和Realizablek-ε模型，并用PIV、micro-PIV等可視化測定結果進行驗證，提出抗堵優化改進模型，魏正英等[12-15,17-23]利用離散相模型模擬離散相致堵顆粒在流道內水沙速度分布、流線分布及沙粒的運動軌跡和濃度的瞬時分布等，但未對不同顆粒的運動特性進行分析研究。

已發表的文獻[15]基于非結構化網格和標準k-e湍流模型開展了迷宮流道內流動和顆粒運動研究，并未對離散相顆粒的運動規律作深入研究；為了深入分析不同性質致堵顆粒在流場內的運動規律，改善迷宮流道滴頭的水力性能，本研究在文獻[15]的基礎上，采用高質量的結構化網格，進行了系列k-e湍流模型的適用性分析，并進一步深入分析了顆粒地運動規律。與已發表的文獻[15]相比，本文更全面考慮顆粒的各種受力，預測不同密度和不同顆粒粒徑下的顆粒運動軌跡變化規律，分析不同顆粒粒徑下顆粒與其周圍流體間的滑移速度變化規律，以及顆粒受到的曳力變化規律等，進一步為灌水器抗堵塞的優化設計提供可靠經驗。

1 物理模型及網格劃分

本文以齒形迷宮流道滴頭為研究對象，進、出口平直流道各延長5 mm，以便于流體的充分發展，取8個流道單元以便于與文獻[15]中的試驗數據進行對照，滴頭結構及參數示意圖如圖1所示。

圖1 滴頭結構及參數示意圖（單位：mm）Fig.1 Structure and parameter of emitter (unit: mm)

應用UG軟件構建物理模型，利用ICEM軟件進行網格劃分，基于分段劃分的思想，將滴頭分為進口段、迷宮流道段和出口段進行結構網格劃分，考慮到流道單元拐角處為三角形塊，選用Y-Block對其進行網格劃分，流道的網格劃分及網格無關性驗證計算如圖 2所示。為了更準確地模擬近壁區的流動狀況，在流道近壁區網格進行加密，并通過網格無關性驗證計算，最終選定網格數量約為27萬個進行數值計算。

圖2 滴頭模型的網格劃分與無關性驗證計算Fig.2 Grid meshing and independent verification calculating of emitter model

2 數值計算與方法分析

為了選擇適宜模擬本文滴頭的湍流模型，開展系列k-e湍流模型的適用性分析，分別采用標準k-ε模型、RNGk-ε模型、Realizablek-ε模型進行滴頭內連續相數值模擬。液相可視為典型的黏性不可壓縮（ρ為常數）流體，并考慮重力作用。

連續相初步收斂后，創建離散相進行固-液相耦合計算。本文采用離散相模型進行稀兩相流數值模擬，在離散相模型中一般要求顆粒相的體積分數小于 10%～12%[24-25]，文中使用單個顆粒，滿足DPM模型使用的要求。

2.1 湍流模型

灌水器流動問題的數值模擬常用的湍流模型有標準k-ε、RNGk-ε和Realizablek-ε模型。為了比較不同湍流模型對迷宮流道滴頭內流動模擬結果的影響，分別采用上述3種k-ε模型進行模擬。

標準k-ε模型是典型的兩方程模型，由湍動能k輸運方程和湍動耗散率ε方程組成。湍動能k是一個關于湍動粘度μt的函數；而RNGk-ε模型與Realizablek-ε模型在標準k-ε模型的基礎上作了進一步改進，通過修正湍動粘度μt將小尺度運動從控制方程中去除，并在湍動耗散率ε方程中添加了反映主流的時均應變率，有效的提高了計算精度，但對于近壁區及低雷諾數流動區域，其模擬存在較大誤差；Realizablek-ε模型進一步在湍動粘度μt的方程上進行改進，引入了有關旋轉和曲率的內容，同時湍動耗散率ε方程中的產生項不再包含k方程中的產生項，有效應用于各種考慮旋轉和大曲率的流動模擬[25-27]。

2.2 離散相模型

離散相模型（discrete phase model，DPM）是兩相流的一種基本模型，其實質是基于歐拉-拉格朗日方法的稀兩相流模型。在求解帶有顆粒的流動中，把流體作為連續相，在歐拉坐標系下處理流場問題，而在拉格朗日坐標系下將顆粒作為離散存在的一個個顆粒。由于實際顆粒運動的隨機性，連續相流場計算初步完成后，需要考慮湍流脈動對顆粒運動的影響，因此本文使用DPM中的隨機軌道模型（stochastic trajectory model）來模擬顆粒物在流場內的運動問題。

對于顆粒在運動過程中受到的作用力，除了考慮顆粒本身的重力外，在顆粒運動過程中顆粒周圍流體在外力的作用下會隨著顆粒作相同方向的運動，相當于顆粒的一部分，稱為虛擬質量力，這部分力對于粒徑不同的顆粒運動有較大的影響。本文綜合考慮重力、浮力、虛擬質量力和曳力的影響，將顆粒的重力和浮力合并表示為FM，對于壓力梯度力、Saffman升力、熱泳力、Basset力等，在本文中忽略不計。力的平衡方程[28-31]如式（1）所示。

式中mp為顆粒質量，kg；u為流體速度，m/s；up為顆粒速度，m/s；Fp為顆粒受到的慣性力，N；FM為考慮了浮力的顆粒相對重力，N；Fd為顆粒受到的曳力，N；Fm為顆粒受到的虛擬質量力，N；Cv為虛擬質量力系數；VP為顆粒體積，m3；dp為顆粒直徑，m；ρp為顆粒密度，kg/m3；g為重力加速度，m/s2；ρ為流體密度，kg/m3；μ為流體動力粘度，Pa·s，取 10–3Pa·s；Rep為顆粒雷諾數（相對雷諾數）；CD為曳力系數，其中系數a1、a2和a3均為常數，按照式（7）計算，式中的系數采用經驗系數[4]，具體取值如表1所示，對于Rep超過1 000的曳力系數，CD近似取為0.44。

表1 曳力系數a1、a2、a3的取值Table 1 constants for drag coefficient a1, a2, a3

2.3 計算方法及模擬參數設置

數值計算采用有限體積法離散控制方程，壓力項等參數采用二階迎風格式，速度和壓力的耦合采用SIMPLE算法求解，流道壁面的處理采用標準壁面函數法，湍流對顆粒的隨機影響采用隨機軌道模型，收斂精度設為10-4。

為準確模擬流道內的顆粒運動，進一步與文獻[15]試驗數據進行對比，模擬中所用的顆粒密度、大小與文獻[15]相同。試驗采用密度分別為1 740 kg/m3的鎂粉粒子、2 500 kg/m3的河砂粒子和2 870 kg/m3的鋁粉粒子以及經過篩網挑選出的直徑在58～65、90～100和125～150μm范圍內的河砂粒子，經過充分攪拌形成濃度為12 g/L的稀兩相流[18]，文獻[15]采用顆粒的直徑上限來開展研究；本文在模擬時，也采用顆粒直徑上限，分別為 65、100和150μm，在進口斷面中心位置水平投入單個顆粒。

3 結果與分析

3.1 湍流模型確定

試驗用水為自來水，試驗壓力為0.01～0.08 MPa，每次測試時間為10 min，測試2次，取平均值計算各灌水器的流量[15]。本文設置進、出口為壓力邊界條件，進口壓力水頭分別為1、2、3、4、5、6、7、8 m，出口壓力水頭為0。

分別利用標準k-ε、RNGk-ε、Realizablek-ε模型對連續相流體進行數值模擬，得到相應的水力性能曲線，如圖3所示。

圖3 水力性能關系曲線Fig.3 Relationship curves of hydraulic performance

滴頭的水力性能利用式（9）來表征，通過式（10）和式（11）計算得到不同模型下水力性能參數k與x，如表2所示。

式中q為出口流量，L/h；h為滴頭的壓力水頭，m；k為滴頭的流量系數；x為滴頭的流態指數。

表2 不同計算模型的水力性能參數k、xTable 2 Hydraulic performance parameters k and x of different calculation models

根據表2中不同工況下的k、x、αk及αx發現，RNGk-ε與Realizablek-ε模型計算得到的k普遍比試驗數據低，x普遍比試驗數據高，標準k-ε模型的k與試驗數據差距很小，相對誤差僅為2.48%。同時，通過對圖3中不同模型的水力性能曲線進行誤差分析發現，RNGk-ε與Realizablek-ε模型模擬得到的水力性能曲線平均模擬誤差分別為8.96%和8.11%，而標準k-ε模型模擬的平均模擬誤差僅為2.32%。

綜合比較發現，標準k-ε模型模擬預測的水力性能關系曲線與文獻[15]試驗數據更為接近，因此后續的離散相顆粒運動模擬基于標準k-ε模型內流場數值計算結果上進行，下文涉及到的數值模擬分析都是基于進口壓力為4 m水頭下進行的。本文的標準k-ε模型模擬誤差略小于文獻[15]的模擬誤差 2.53%，這主要是由于文獻[15]采用非結構化網格且網格數量較少，而本文采用高質量的結構化網格，為了保證后續離散相模擬過程的準確性，對計算模型進行高質量的網格劃分，能進一步減少模擬誤差，使模擬預測結果更為接近實際。

3.2 連續相分析

圖4為流道縱向對稱軸面的流場壓力云圖，觀察圖4發現，流道內的壓力分布沿水流運動方向逐漸下降且變化均勻，壓力變化基本發生在流道齒尖處，水流流經每個流道單元壓力下降值基本相等。

圖4 流道縱向對稱軸面流場壓力云圖Fig.4 Pressure contour of longitudinal symmetry plane

圖 5為流道不同截面的速度云圖、等值線圖和流線圖。觀察圖5a發現，速度在每個流道單元內的分布規律相同；觀察圖5b中第5流道單元的速度云圖和流線圖發現，流道內流場按照流速大小可以分為 3部分，分別為低速區A、主流區B和高速發展區C，低速區A位于背水面齒尖拐角處，流線呈漩渦狀，也稱漩渦區[12]，主流區B流線平滑密集流速較高；圖5c為1-1橫截面速度等值線圖，速度等值線關于y軸對稱分布，從流道邊壁到齒尖位置速度等值線圖越來越密集，說明越靠近齒尖速度變化越劇烈，能量耗散越多。

圖5 流道不同截面的速度云圖、等值線圖Fig.5 Velocity nephogram and velocity isolines contour of different sections in flow channel

綜合圖5b和5c發現，高速發展區C靠近齒尖位置，流速在0.2～2.6 m/s范圍變化，低速區A位于流道邊壁附近，流速在0.2～0.8 m/s范圍內變化，主流區B流速在0.8～2.0 m/s范圍內變化。從而得出，流場內能量消耗主要發生在來流側齒尖處，流體的動能轉化為其他形式的能量耗散掉，從而達到消能效果。因此，為了提高滴頭的消能能力，通過研究齒尖處拐角大小和結構形式與其他結構參數之間的關系來進一步設計出消能效果最好的流道結構。

3.3 離散相分析

3.3.1 不同密度顆粒的運動分析

分別模擬直徑同為100μm，密度為1 740 kg/m3的鎂粉、2 500 kg/m3的河砂和2 870 kg/m3的鋁粉單個粒子在流場內的運動，結果如圖6所示。

對比圖6a、6d和6g中的流場速度云圖可以看出，隨著顆粒密度的增大，轉角附近的高速發展區逐漸消失，主流區與低速區開始互相摻混，流場對顆粒的擾動作用增大。

對比圖6b、6e和6h中顆粒模擬運動軌跡發現，顆粒直徑相同的情況下，隨著顆粒密度的增大，運動軌跡在旋渦區內路程增多。鎂粉、河砂和鋁粉粒子通過齒形迷宮流道滴頭經過的路程分別為43.97，44.25和50.79 mm，由此可知，隨著顆粒密度的增大，顆粒通過流道的路程增多，顆粒運動跟隨性變差。圖6b、6e和6h中粒子流經1-1橫截面時，所對應的速度分別為1.09、1.05和1.00 m/s，由此可知，顆粒運動跟隨性變差會導致顆粒速度減小。

對比圖6b、6e和6h中CFD模擬軌跡與圖6c、6f和6i中文獻[15]PIV試驗觀測到的運動軌跡發現，CFD模擬得到的顆粒運動軌跡與 PIV試驗觀測到的顆粒運動軌跡雖然有所不同，但是軌跡運動的變化規律相同，這是由于文獻[15]PIV試驗測量結果的顆粒直徑和顆粒釋放的準確位置無法確定，故與本研究的模擬結果存在較小差異。

3.3.2 不同直徑顆粒的運動分析

分別模擬顆粒密度同為2 500 kg/m3，直徑為65、100和150μm的單個河砂粒子在流道內的運動，得到結果如圖7所示。

對比圖7a、7b和7c中顆粒的運動軌跡發現，當顆粒直徑為65μm時，顆粒隨主流區的水流運動，運動軌跡光滑有規律；當直徑為100μm時，顆粒與流道邊壁發生多次碰撞，顆粒的運動軌跡開始不再光滑，但運動軌跡仍然沿著主流區；當直徑為150μm時，顆粒脫離主流區不斷進入漩渦區，出現小漩渦現象，運動軌跡沒有規律。3種直徑顆粒通過流道的路程分別為 34.04，44.25和49.06 mm，由此可知，隨著顆粒直徑的增大，顆粒在漩渦區運動的路程增多，顆粒運動跟隨性變差，顆粒直徑變化較密度變化對顆粒運動影響大。

對比圖7d、7e和7f中砂粒與周圍流體的速度變化發現，當直徑為65μm時，顆粒與流體的速度差異不大，顆粒與水流運動保持良好的一致性，速度整體變化幅值為3.60 m/s，速度均值為1.42 m/s；當直徑為100μm時，顆粒與流體的速度差異增大，顆粒運動與流體基本保持一致，速度整體變化幅值為 2.19 m/s，速度均值為1.14 m/s；當直徑為150μm時，顆粒速度與流體速度偏離更為明顯，并出現了固液速度錯峰的現象，速度變化越來越沒有規律，速度整體變化幅值為1.93 m/s，速度均值為1.10 m/s。由此得出，隨著顆粒直徑的增大，顆粒與周圍流體之間的速度偏差不斷增大，速度整體變化幅值與均值減小。這是因為隨著顆粒直徑增大，顆粒運動過程中受到的阻力增大，使得顆粒運動跟隨性變差，速度減小。

圖6 不同密度顆粒所在流場的速度云圖及其運動軌跡Fig.6 Flow field velocity nephogram and moving trajectory of single particle with different densities

圖7 不同直徑砂粒的運動軌跡、與周圍流體的速度變化對比圖Fig.7 Moving trajectory and velocity change with surrounding fluid of single sand particle with different diameters

圖 8為不同直徑砂粒所受曳力對比圖，顆粒受到的曳力越來越大，變化幅度明顯，結合圖7a、7b和7c中顆粒的運動軌跡發現，隨著顆粒直徑增大，顆粒受到的曳力增大，顆粒運動隨流性變差。

圖8 不同直徑單個砂粒所受曳力對比Fig.8 Drag force contrast of single sand particle with different diameters

3.3.3 顆粒運動速度變化分析

圖9和圖10分別為在不同密度和不同直徑下，本文與文獻[15]對相同的灌水器模型進行 CFD模擬和相應試驗結果的比較，橫坐標代表顆粒在流道長度方向的位置。圖9中，隨著顆粒密度增大，在13.0～16.9 mm這段流道長度內顆粒達到極值的次數分別為8、7和6次；圖10中，隨著顆粒直徑增大，在13.0～17.4 mm這段流道長度內顆粒達到極值的次數分別為9、7和6次。相同長度的流道單元內速度極值出現的次數減少，顆粒速度變化變劇烈，這是由于隨著密度或直徑的增大，顆粒受到阻力的影響增大，在慣性的作用下易被甩離高速發展區，速度減小，顆粒逐漸往速度減小的區域運動。

圖11為不同工況下的顆粒速度幅值變化圖，為了保證結果的準確度，對圖11中每個工況下的速度都取平均值。圖11a和11b分別為顆粒運動速度極大值隨著顆粒密度和直徑增大的變化圖，速度極大值隨著密度或直徑的增大變化不大，均在2.1～2.5 m/s內變化，本文運動速度極大值的模擬誤差為 1.5%，而文獻[15]的模擬誤差為6.0%；圖11c和11d分別為顆粒運動平均速度隨著顆粒密度和直徑增大的變化圖，平均速度隨著密度或直徑的增大逐漸減小，在1.1～1.8 m/s內變化，減小幅度基本相同，本文運動平均速度的模擬誤差為 2.7%，而文獻[15]的模擬誤差為7.1%；圖11e和11f分別為顆粒運動速度極小值隨著顆粒密度和直徑增大的變化圖，隨著顆粒直徑和密度增大，速度幅值減小程度從大到小依次為：顆粒直徑、顆粒密度，本文運動速度極小值的模擬誤差為3.0%，而文獻[15]的模擬誤差為 9.2%。由此說明，隨著顆粒直徑或密度增大，顆粒整體運動速度減小，顆粒運動速度極小值受其影響最大，顆粒運動平均速度受其影響次之，顆粒運動速度極大值受其影響最小。

圖9 不同密度的顆粒運動速度對比Fig.9 Comparison of particle velocity with different densities

圖10 不同直徑的砂粒運動速度對比Fig.10 Comparison of sand particle velocity with different diameters

圖11 不同工況下的顆粒運動速度幅值對比Fig.11 Comparison of particle movement velocity amplitude in different conditions

圖 11中本文對顆粒運動速度的平均模擬誤差為2.34%，比文獻[15]平均模擬誤差減少了 5.05%，結果更加接近 PIV試驗數據，這是由于本文針對顆粒受力不僅考慮了重力與曳力，還考慮了顆粒的浮力與虛擬質量力對于離散顆粒運動的影響，比文獻[15]只考慮重力與曳力的模擬結果更符合顆粒運動客觀規律。

4 結 論

1）通過k-ε系列模型的數值模擬結果與PIV試驗數據的對比，分析了各湍流模型在迷宮流道滴頭內流動的適應性，結果表明，標準k-ε模型預測結果最接近PIV試驗數據，水力性能曲線平均模擬誤差僅為2.32%。

2）隨著顆粒密度增大，液相流場受顆粒擾動變大且主流區與低速區開始互相摻混，顆粒運動跟隨性變差；隨著顆粒直徑的增大，顆粒與周圍流體之間的速度偏差不斷增大，速度整體變化幅值與均值減小。顆粒受到的曳力增大，顆粒運動速度與液相速度出現峰值錯亂現象，顆粒運動軌跡變亂，顆粒在慣性力的作用下更容易被帶入漩渦區，顆粒的運動速度減小，最后可能沉積在流道拐角內側的漩渦區，易造成迷宮流道堵塞。

3）隨著顆粒直徑或密度增大，顆粒整體運動速度減小，顆粒運動速度極小值受影響最大，顆粒平均運動速度受影響次之，顆粒運動速度極大值受影響最小，從而證明，級配顆粒運動速度隨其直徑和密度不同而變化很大。

4）通過隨機軌道模型模擬顆粒運動，基于文獻[15]的數值研究，進一步綜合考慮了虛擬質量力、重力、曳力和浮力對于顆粒運動的影響，模擬結果更加接近 PIV試驗數據，顆粒運動速度的平均模擬誤差僅為2.34%，顆粒運動軌跡更符合客觀規律。

[1]張俊，趙萬華，粟曉玲，等. 微灌長流道灌水器結構特性的研究綜述[J]. 農業工程學報，2005，21(1)：182－185.Zhang Jun, Zhao Wanhua, Su Xiaoling, et al. Review of structural characteristics analysis of the long-path emitters formicroirrigation[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2005,21(1): 182－185. (in Chinese with English abstract)

[2]Sahin U, Anapali O, D?nmez M F, et al. Biological treatment of clogged emitters in a drip irrigation system[J]. Journal of Environmental Management, 2005, 76(4): 338－341.

[3]Wei Q, Gang L, Jie L, et al. Evaluations of emitter clogging in drip irrigation by two-phase flow simulations and laboratory experiments[J]. Computers & Electronics in Agriculture, 2008, 63(2): 294－303.

[4]Zhang Jun, Zhao Wanhua, Tang Yiping, et al. Anti-clogging performance evaluation and parameterized design of emitters with labyrinth channels[J]. Computers and Electronics in Agriculture, 2010, 74(1): 59－65.

[5]馮吉，孫昊蘇，李云開. 滴灌灌水器內顆粒物運動特性的數字粒子圖像測速[J]. 農業工程學報，2013，29(13)：90－97.Feng Ji, Sun Haosu, Li Yunkai. Visualizing particles movement characteristics in drip irrigation emitters with digital particle image velocimetry[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2013, 29(13): 90－97. (in Chinese with English abstract).

[6]李云開，馮吉，宋鵬，等. 低碳環保型滴灌技術體系構建與研究現狀分析[J]. 農業機械學報，2016，47(6)：83－92.Li Yunkai, Feng Ji, Song Peng, et al. Developing situation and system construction of low-carbon environment friendly drip irrigation technology[J]. Transactions of the Chinese Society for Agricultural Machinery, 2016, 47(6): 83－92. (in Chinese with English abstract)

[7]李云開，馮吉. 滴灌灌水器內部水動力學特性測試研究進展[J]. 排灌機械工程學報，2014，32(1)：86－92.Li Yunkai, Feng Ji. Progress in measurement of hydrodynamic characteristics in drip irrigation emitters[J].Journal of Drainage and Irrigation Machinery Engineering,2014, 32(1): 86－92. (in Chinese with English abstract)

[8]李云開，楊培嶺，任樹梅，等. 圓柱型灌水器迷宮式流道內部流體流動分析與數值模擬[J]. 水動力學研究與進展(A 輯)，2005，20(6)：736－743.Li Yunkai, Yang Peiling, Ren Shunmei, et al. Analyzing and modeling flow regime in labyrinth path drip irrigation column emitter with CFD[J]. Journal of Hydrodynamics: Ser A, 2005, 20(6): 736－743. (in Chinese with English abstract)

[9]李云開，楊培嶺，任樹梅，等. 圓柱型迷宮式流道滴灌灌水器平面模型試驗研究[J]. 農業機械學報，2006，37(4)：48－51.Li Yunkai, Yang Peiling, Ren Shunmei, et al. Plane model of labyrinth path drip irrigation column emitter[J]. Transactions of the Chinese Society for Agricultural Machinery, 2006,37(4): 48－51. (in Chinese with English abstract)

[10]Liu H, Huang G. Laboratory experiment on drip emitter clogging with fresh water and treated sewage effluent. Agric Water Manag[J]. Agricultural Water Management, 2009,96(5): 745－756.

[11]魏正英，趙萬華，唐一平，等. 滴灌灌水器迷宮流道主航道抗堵設計方法研究[J]. 農業工程學報，2005，21(6)：1－7.Wei Zhengying, Zhao Wanhua, Tang Yiping, et al.Anti-clogging design method for the labyrinth channels of drip irrigation emitters[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2005, 21(6): 1－7. (in Chinese with English abstract)

[12]魏正英，唐一平，溫聚英，等. 灌水器微細流道水沙兩相流分析和微PIV及抗堵實驗研究[J]. 農業工程學報，2008，24(6)：1－9.Wei Zhengying, Tang Yiping, Wen Juying, et al. Two-phase flow analysis and experimental investigation of micro-PIV and anti-clogging for micro-channels of emitter[J].Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2008, 24(6): 1－9.(in Chinese with English abstract)

[13]Wei Zhengying, Cao Meng, Liu Xia, et al. Flow behaviour analysis and experimental investigation for emitter micro-channels[J]. Chinese Journal of Mechanical Engineering, 2012, 25(4): 729－737.

[14]喻黎明，梅其勇. 迷宮流道灌水器抗堵塞設計與PIV試驗[J]. 農業機械學報，2014，45(9)：155－160.Yu Liming, Mei Qiyong. Anti-clogging design and experimental investigation of PIV for labyrinth-channel emitters of drip irrigation emitters[J]. Transactions of the Chinese Society for Agricultural Machinery, 2014, 45(9): 155－160. (in Chinese with English abstract)

[15]喻黎明，吳普特，牛文全，等. 迷宮流道內固體顆粒運動的CFD模擬及PIV驗證[J]. 農業機械學報，2009，40(5)：45－51.Yu Liming, Wu Pute, Niu Wenquan, et al. CFD numerieal simulation and PIV verification about the movement of solid particles in labyrinth channel[J]. Transactions of the Chinese Society for Agricultural Machinery, 2009, 40(5): 45－51. (in Chinese with English abstract)

[16]牛文全，喻黎明，吳普特，等. 迷宮流道轉角對灌水器抗堵塞性能的影響[J]. 農業機械學報，2009，40(9)：51－55.Niu Wenquan, Yu Liming, Wu Pute, et al. Influence of angle of labyrinth channels on anti-clogging performance of emitter[J]. Transactions of the Chinese Society for Agricultural Machinery, 2009, 40(9): 51－55. (in Chinese with English abstract)

[17]張俊, 魏公際, 趙萬華. 灌水器內圓弧形流道的液固兩相流場分析[J]. 中國機械工程, 2007, 18(5)：589－593.Zhang Jun, Wei Gongji, Zhao Wanhua. Numerical analysis on liquid-solid two-phase flows in arc-type channel of the emitter[J]. China Mechanical Engineering, 2007, 18(5): 589－593. (in Chinese with English abstract)

[18]喻黎明. 滴灌灌水器水力特性及抗堵塞研究[D]. 楊凌：西北農林科技大學，2009.Yu Liming. Study on Hydraulic and Anti-clogging Performance of Drip Emitter[D]. Yangling: Northwest A&F University, 2009. (in Chinese with English abstract)

[19]常瑩華，牛文全，王維娟. 滴灌灌水器迷宮流道的內部流體數值模擬與流動分析[J]. 西北農林科技大學學報(自然科學版)，2009，37(2)：203－208.Chang Yinghua, Niu Wenquan, Wang Weijuan. Numerical simulation and flow analysis of labyrinth path of drip irrigation emitters[J]. Journal of Northwest A&F University(Natural Science Edition), 2009, 37(2): 203－208. (in Chinese with English abstract)

[20]喻黎明，譚弘，鄒小艷，等. 基于CFD-DEM耦合的迷宮流道水沙運動數值模擬[J]. 農業機械學報，2016，47(8)：65－71.Yu Liming, Tan Hong, Zou Xiaoyan, et al. Numerical simulation of water and sediment flow in labyrinth channel based on coupled CFD-DEM[J]. Transactions of the Chinese Society for Agricultural Machinery, 2016, 47(8): 65－71. (in Chinese with English abstract)

[21]鄧濤，魏正英，王立朋，等. 壓力補償灌水器分步式CFD設計與實驗研究[J]. 清華大學學報(自然科學版)，2012，52(4)：513－516.Deng tao, Wei Zhengying, Wang lipeng, et al. Stepwise CFD simulation and experiments of pressure-compensating emitters[J]. Journal of Tsinghua University (Natural Science Edition), 2012, 52(4): 513－516. (in Chinese with English abstract)

[22]金文，張鴻雁，何文博. 灌水器流道結構對水力性能影響的數值分析[J]. 中國農業大學學報，2012，17(2)：139－143.Jin Wen, Zhang Hongyan, He Wenbo, et al. Numerical investigation on effect of channel structure to hydraulic performance of emitter[J]. Journal of China Agricultural University, 2012, 17(2): 139－143. (in Chinese with English abstract)

[23]張重，任樹梅，楊培嶺，等. 圓柱型灌水器迷宮式流道內部流體流動的CFD模擬精度研究[J]. 中國農業大學學報，2016，21(9)：149－155.Zhang Chong, Ren Shumei, Yang Peiling, et al. Research on the CFD simulation accuracy of labyrinth path drip irragation cylindrical emitter[J]. Journal of China Agricultural University, 2016, 21(9): 149－155. (in Chinese with English abstract)

[24]李亞林，袁壽其，湯躍，等. 離心泵內示蹤粒子運動的離散相模型模擬[J]. 農業機械學報，2012，43(11)：113－118.Li Yalin, Yuan Shouqi, Tang Yue, et al. Simulation of tracer particles movement by discrete phase model in the centrifugal pump[J]. Transactions of the Chinese Society for Agricultural Machinery, 2012, 43(11): 113－118. (in Chinese with English abstract)

[25]唐學林，李林偉，蘇硯文，等. 螺旋槽流道微泵的數值模擬方法分析[J]. 排灌機械工程學報，2010，28(1)：31－37.Tang Xuelin, Li Linwei, Su Yanwen, et al. Analysis of numerical simulation method for flow inside micro-machined spiral-channel pump[J]. Journal of Drainage and Irrigation Machinery Engineering, 2010, 28(1): 31－37. (in Chinese with English abstract)

[26]金文. 微灌滴頭微通道內流流場實驗研究與數值模擬[D].西安：西安建筑科技大學，2010.Jin Wen. Experimental Research and Numerical Simulating Approaches on Internal Flow Fields of Micro Channel in Micro-irrigation Emitter[D]. Xi￠an: Xi'an University of Architecture and Technology, 2010. (in Chinese with English abstract)

[27]王福軍. 計算流體動力學分析：CFD軟件原理與應用[M].北京：清華大學出版社，2004.

[28]Zhang J, Zhao W, Tang Y, et al. Numerical investigation of the clogging mechanism in labyrinth channel of the emitter[J].International Journal for Numerical Methods in Engineering,2010, 70(13): 1598－1612.

[29]Fluent Inc. Fluent6. 3 User's Guide[M]. Lcbanon: Fluent Inc.,2006.

[30]Mei R, Adrian R J, Hanratty T J. Particle dispersion in isotropic turbulence under Stokes drag and Basset force with gravitational settling[J]. Journal of Fluid Mechanics, 1991,225: 481－495.

[31]陳彬，顏歡，劉閣，等. 流固兩相流的稀疏離散相模型研究進展[J]. 化工進展，2016，35(11)：3400－3412.Chen Bin, Yan Huan, Liu Ge, et al. Progress in the dilute discrete model of fluid-solid two phase flow[J]. Chemical Industry and Engineering Progress, 2016, 35(11): 3400－3412.(in Chinese with English abstract)