導流板對青花椒堆積式烘干烘房內流場均勻性的影響

譚均， 王教領， 王子軒， 劉恩谷，黃園園， 楊明金， 楊玲， 蒲應俊

1. 丘陵山區農業裝備重慶市重點實驗室/西南大學 工程技術學院，重慶 400715；2. 農業農村部南京農業機械化研究所，南京 210014

青花椒Zanthoxylumschinifolium為蕓香科植物, 是十分重要的經濟作物[1]． 青花椒是我國調味調料品的重要來源, 且有藥用價值, 廣泛種植于四川、 重慶、 云南、 貴州等地[2]． 近年來, 隨著農業產業結構調整, 其種植面積不斷擴大[3]． 新鮮青花椒含水率高, 在儲存和運輸過程中容易因霉爛、 褐變而使品質降低． 對青花椒進行烘干處理可降低其含水率、 延長儲存時間、 提升經濟價值[4]． 傳統堆積式烘干采用挖坑并結合磚砌構建烘房, 坑部形成高度為400～500 mm的勻風室, 上部堆積厚度為1.2～1.5 cm帶枝花椒或凈椒, 烘房上部敞開或以麻袋覆蓋, 熱源多為燃煤式熱風爐, 能耗高、 污染重、 安全性差, 屬于需淘汰的落后生產工藝裝備． 采用空氣能熱泵作為替代熱源, 用工業絕熱板材制成閉式烘房, 自動控制水分排出, 具有能效高、 環境友好、 安全性好、 控制方便等優點, 逐漸成為青花椒烘干領域的研究熱點[5]． 烘房內流場均勻性決定了烘干速率、 烘干品質和能耗． 國內外研究人員借助CFD(Computational Fluid Dynamics)對烘房內流場均勻性進行了相關研究, 主要集中在改進烘房結構、 增設輔助勻風裝置兩方面． 陳忠加等[6]、 牛天寶等[7]探究了不同送風回風方式對烘房內部溫度、 速度均勻性的影響; 師建芳等[8]、 張健平[9]探究了烘房入口處不同風機數量對烘房內流場均勻性的影響; 龔中良等[10]、 代建武等[11]通過改變出風孔徑、 排列方式等來優化其內部流場; 劉瑞等[12]針對氣體射流沖擊干燥機氣流分配室提出了3種不同的優化結構; Wang等[13]通過在烘房進風口處安裝引風罩優化烘房內流場分布; Song等[14]研究了入口氣流速度、 入口直徑、 隔板數量對西洋參干燥室內氣流和溫度分布的影響． 學者們還探究了烘房內不同形狀、 不同數量的導流板對其內部流場均勻性的影響[15-20]; 王鵬浩[21]探究了排風速度、 導流板角度、 導流板位置對空載和加載工況烘房內流場的影響; 陳竹筠[22]探究了不同間距的格柵對熱泵型海帶烘房內氣流的影響． 綜上可知, 烘房結構及導流板的設置都會影響其流場均勻性．

本文針對青花椒空氣能熱泵烘干裝備, 采用數值研究和實驗研究相結合的手段, 探討導流板對青花椒堆積式烘干烘房內流場均勻性的影響, 基于BBD試驗和統計分析優化了導流板尺寸和工作參數, 可為青花椒高品質綠色烘干工藝及烘房設計提供參考．

1 材料與方法

1.1 材料與儀器設備

材料: 青花椒(九葉青), 2022年6月購于重慶市璧山區; 細食用鹽, 用于密度法測量孔隙率．

儀器設備: 5HGKB-4型空氣能熱泵烘干機(自制); WD412型風速傳感器(杭州且遠科技有限公司); FBS-760A型鹵素水分測定儀(廈門市弗布斯檢測設備有限公司); BSA224S-CW型電子天平(賽多利斯科學儀器(北京)有限公司); 燒杯、 量筒若干．

1.2 方法

1.2.1 孔隙率測量和烘房加載

新鮮青花椒含水率為67.54%, 鮮果直徑5.5 mm． 參考文獻[23]采用密度法測量孔隙率, 如圖1所示． 青花椒堆積密度和真實密度分別為641.79 kg/m3和961.59 kg/m3, 計算得到青花椒孔隙率為33.0%． 烘房加載方式如圖2所示, 裝載厚度為30 mm．

圖1 孔隙率測量

1.2.2 仿真模型建立

空氣能熱泵烘干機外形如圖3所示, 其烘房由出風口、 烘干區、 門、 勻風室、 導流板、 勻風板、 進風口等組成, 如圖4所示． 采用SolidWorks建立烘房三維模型, 主要尺寸參數為: 內腔長寬高1 800 mm×1 300 mm×1 800 mm; 進風口330 mm×220 mm; 出風口700 mm×200 mm; 勻風板1 800 mm×1 300 mm; 勻風板孔徑3 mm、 厚度3 mm、 孔隙率18.1%; 勻風室高度330 mm; 門尺寸1 100 mm×1 500 mm(模型中當作壁面)． 導流板尺寸根據實驗方案確定．

圖4 烘房結構

建模時, 為便于數值模擬和實驗驗證, 進行如下處理: ① 不考慮烘房壁厚, 絕熱壁面; ② 不考慮勻風板厚度, 并設置無厚度多孔躍升邊界條件[24]; ③ 干燥介質空氣視為不可壓縮理想氣體[25]．

采用無厚度多孔躍升邊界條件, 可避免考慮勻風板上約6萬個孔徑為3 mm小孔帶來的計算資源消耗問題． 30 mm裝載厚度屬于深床干燥范疇, 可考察堆積式烘干在一定裝載量條件下的烘房流場均勻性．

采用ICEM(International Council for Educational Media)軟件對烘房三維模型進行非結構網格劃分, 如圖5所示． 為保證計算精度, 進風口、 出風口、 勻風板及物料層的網格密度大于壁面． 空載烘房網格數量為2.67×106, 加載烘房網格數量為2.8×106． 在進行網格無關性和模型可靠性驗證時未考慮導流板; 在進行單因素試驗和BBD試驗時, 考慮導流板(圖6)． 試驗導流板參數包括導流板角度θ、 導流板位置S、 導流板長度L． 固定導流板其他參數, 包括寬度300 mm、 孔徑15 mm、 孔間距60 mm、 厚度1 mm．

圖5 烘房三維模型網格劃分

圖6 導流板安裝及結構

烘房邊界條件及參數設置見表1． 進風口選擇速度入口, 其大小由風速傳感器測量得到; 出風口選擇壓力出口; 參考文獻[27]計算入口、 出口湍流強度及水力直徑; 勻風板選擇多孔躍升[25]; 物料層設置為多孔介質, 參考文獻[24]方法測量孔隙率、 黏性阻力系數及慣性阻力系數; 烘房其他部分選擇無滑移絕熱壁面邊界條件．

表1 邊界條件及參數

采用Ansys Fluent軟件對三維模型求解, 選擇Realizable k-ε湍流模型及SIMPLE算法． 壓強設置為標準大氣壓, 在豎直方向重力加速度設置為-9.81 m/s2, 采用混合初始化方式, 收斂精度設置為1×10-4． Realizable k-ε模型中湍動能k和擴散率ε的方程分別為:

(1)

(2)

式(1)、 式(2)中,ρ為密度, kg/m3;μ為層流黏性系數, Pa·s;μt為湍流黏性系數, Pa·s;xi,xj分別為位移在i,j方向上的分量, m;Gk為層流速度梯度產生的湍流動能, m2/s2;Gb為浮力產生的湍流動能, m2/s2;C1g,C2g,C3g為湍流經驗常數;σk為湍動能的湍流普朗特常數;σε為耗散率的湍流普朗特常數． 模型參數設置同文獻[22]．

1.2.3 試驗方法

通過單因素試驗分析導流板角度、 導流板位置和導流板長度對烘房內流場及流場均勻性的影響, 并確定BBD試驗因子水平． 通過BBD試驗確定各試驗因子的影響程度, 并優化導流板尺寸和工作參數． 表2為單因素試驗安排, 試驗號5,8,15具有相同的因子組合, 僅試驗1組． 表3為BBD試驗因子水平．

1.2.4 評價指標

根據速度云圖和速度不均勻系數評價烘房流場均勻性．

速度不均勻系數為監測平面上觀測點速度均方差與平均速度的比值[26], 可表達為:

(3)

(4)

實測時, 選取3個監測平面, 其z坐標分別為600 mm、 1 000 mm、 1 400 mm, 每個平面選取9個觀測點, 選取位置從平面中心向x和y方向雙向展開,x方向間距為400 mm,y方向間距為250 mm, 即x方向布置3行,y方向布置3列．x,y,z坐標設置見圖6． 測量時, 每個觀測點測量3組數據, 每組數據100次重復． 仿真時, 選取5個監測平面, 其z分別為600 mm、 800 mm、 1 000 mm、 1 200 mm、 1 400 mm, 每個平面選取377個觀測點, 選取位置從平面中心向x和y方向雙向展開,x方向間距為60 mm,y方向間距為100 mm, 即x方向布置29行,y方向布置13列．

2 結果與分析

2.1 網格無關性及模型可靠性驗證

選擇不同網格劃分策略改變網格大小, 得到1 598 302,2 369 430,2 665 263,4 001 428共4種網格數目進行網格無關性驗證． 空載時在烘房坐標系中y和z分別為1 000 mm的兩個平面交線上, 驗證不同網格劃分策略下風速隨x坐標位置的變化情況, 如圖7所示． 由圖7可知, 不同的網格劃分策略風速隨x坐標位置變化趨勢基本相同, 當網格數目為2 665 263時, 增加網格數目對風速變化影響很?。?因此, 為了減少計算量, 采用網格數目為2 665 263的網格劃分策略進行后續數值研究．

圖7 網格數目對風速的影響

分別選取z坐標為600 mm、 1 000 mm、 1 400 mm 3個平面上的風速進行模型可靠性驗證． 空載時3個平面各觀測點風速模擬值和實驗值數據見表4． 3個平面風速模擬值和實驗值相對誤差的平均值分別為5.87%,5.78%,7.38%． 根據各觀測點風速數據, 可以計算出3個平面上基于模擬風速和實測風速的速度不均勻系數, 如圖8a所示． 3個平面上模擬風速的速度不均勻系數分別為0.48,0.61,0.58, 實測風速的速度不均勻系數分別為0.46,0.6,0.61, 其相對誤差分別為4.35%,1.67%,4.92%, 平均為3.65%． 加載時, 3個平面各觀測點風速模擬值和實驗值數據見表5． 3個平面風速模擬值和實驗值相對誤差的平均值分別為9.87%,9.24%,9.16%． 3個平面上基于模擬風速和實測風速的速度不均勻系數, 如圖8b所示． 3個平面上模擬風速的速度不均勻系數分別為0.42,0.87,0.67, 實測風速的速度不均勻系數分別為0.41,0.85,0.68, 其相對誤差分別為2.44%,2.35%,1.47%, 平均為2.09%． 風速和速度不均勻系數的試驗和仿真結果表明, 本文建立的仿真模型可靠性高, 可用于后續相關研究．

表4 模型可靠性驗證(空載)

表5 模型可靠性驗證(加載)

圖8 監測平面速度不均勻系數

2.2 單因素試驗

圖9和圖10分別為單因素試驗風速的速度云圖和速度不均勻系數影響圖．

圖9 單因素試驗速度云圖(y=650 mm豎直面)

隨著導流板角度增大, 導流板對進風口高速流動氣流的擋風效果增強, 一部分氣流流向烘干區, 能夠很大程度地避免勻風室內高速流體直接流向右側(門側), 然后經上部空間直接回流出風口(圖9a)． 與導流板角度30°,45°相比, 導流板角度為60°,75°,90°時的速度不均勻系數更低(圖10a)． 因此, 選定導流板角度為60°,75°,90°進行后續BBD試驗．

隨著導流板位置的改變, 導流板對進風口高速流動氣流沿烘房空間分布作用明顯, 導流板太靠近進風口(如位置500 mm處)或太靠近右側(如位置1 300 mm處)都不利于流場均勻性的提高． 導流板太靠近進風口, 在右側上方區域形成低速區; 導流板太靠近右側, 在烘房中心區域形成低速區(圖9b)． 與導流板位置500 mm、 1 300 mm相比, 導流板位置為700 mm、 900 mm、 1 100 mm時的速度不均勻系數更低(圖10b)． 因此, 選定導流板位置為700 mm、 900 mm、 1 100 mm進行后續BBD試驗．

隨著導流板長度的增加, 導流板對進風口高速流動氣流的擋風效果增強, 使經導流板兩側流向烘房右側的氣流減少, 速度云圖均勻性降低(圖9c)． 與導流板長度1 100 mm、 1 300 mm相比, 導流板長度為500 mm、 700 mm、 900 mm時的速度不均勻系數更低(圖10c)． 但是, 當監測平面高度增大到1 200 mm后, 導流板長度為500 mm時的速度不均勻系數較大, 且監測平面平均速度較低． 綜合考慮, 選定導流板長度為700 mm、 900 mm、 1 100 mm進行后續BBD試驗．

2.3 BBD試驗

使用Design-Expert 12中的Box-Behnken模塊, 以導流板角度、 導流板位置和導流板長度作為試驗因子, 以z=1 000 mm平面上的速度不均勻系數為響應值, 設計3因子3水平BBD試驗． 表6為BBD試驗設計與結果, 表7為BBD試驗方差分析．

表6 BBD試驗設計與結果

表7 BBD試驗方差分析

由統計分析可知, A,AB,A2,B2和C2對速度不均勻系數影響極具有統計學意義, B,C對速度不均勻系數影響具有統計學意義(p<5%), 其他項對速度不均勻系數影響不具有統計學意義． 各單因素對速度不均勻系數影響從大到小依次為導流板角度、 導流板長度、 導流板位置．

對表7中數據進行回歸擬合, 剔除不具有統計學意義項, 得到二次回歸方程:

(5)

圖11為速度不均勻系數交互作用響應面圖． 隨著導流板角度及導流板與進風口距離增大, 速度不均勻系數呈現出先減小后緩慢增加的趨勢, 導流板角度90°、 導流板位置900～1 000 mm區間具有較低的速度不均勻系數(圖11a); 隨著導流板與進風口距離及導流板長度增大, 速度不均勻系數呈現出先緩慢減小的后緩慢增加的趨勢, 導流板長度900 mm、 導流板位置900 mm區域具有較低的速度不均勻系數(圖11b)． 圖11c所示, 較低速度不均勻系數對應導流板角度與導流板長度參數與圖11a和圖11b一致．

圖11 速度不均勻系數交互作用響應面

2.4 烘房優化

根據BBD試驗回歸方程, 以速度不均勻系數最小為目標, 得到最優導流板參數為導流板角度86.01°、 導流板位置924.75 mm、 導流板長度855.51 mm． 此時,z=1 000 mm平面上的速度不均勻系數預測值為0.384 9． 根據導流板最優參數, 建立基于CFD的烘房三維模型并求解, 可以得到不同平面上的速度不均勻系數模擬值． 其中,z=1 000 mm平面上的速度不均勻系數模擬值為0.3887． 預測值與模擬值之間的相對誤差為0.98%． 分別選取z為600 mm、 1 000 mm、 1 400 mm 3個平面上的風速進行仿真與試驗對比． 優化烘房空載時, 3個平面各觀測點風速模擬值和實驗值數據見表8． 3個平面風速模擬值和實驗值相對誤差的平均值分別為6.03%,6.43%,9.08%． 優化烘房加載時, 3個平面各觀測點風速模擬值和實驗值數據見表9． 3個平面風速模擬值和實驗值相對誤差的平均值分別為9.98%,8.70%,6.83%． 因此, 該仿真結果與烘房內部實際結果相近, 可用于探究烘房內部的流場變化．

表8 優化烘房仿真與試驗結果對照(空載)

表9 優化烘房仿真與試驗結果對照(加載)

根據導流板最優參數, 建立烘房三維模型, 按照前述方法求解得到優化烘房流場速度云圖, 并與相同參數原始烘房(未安裝導流板)進行比較, 如圖12所示． 監測優化烘房和原始烘房截面高度600 mm、 800 mm、 1 000 mm、 1 200 mm、 1 400 mm各觀測點風速, 從而得到各水平面上的速度不均勻系數(表10), 并繪出速度不均勻系數隨監測高度的變化曲線, 如圖13所示．

表10 烘房流場均勻性指標

圖12 優化烘房和原始烘房速度云圖

圖13 烘房速度不均勻系數

由圖12可知, 優化烘房流場受到導流板的作用, 其上部烘干區域形成2個子區間, 有利于烘干區域流場均勻性的提高, 烘房加載后其流場均勻性有所降低． 由表10和圖13可知, 與原始烘房相比, 優化烘房(空載)各截面高度平均速度不均勻系數降低39.2%, 優化烘房(加載)各截面高度平均速度不均勻系數降低20.1%, 烘房內流場均勻性得到較大提高．

3 結論

1) 基于Realizable k-ε湍流模型和勻風板多孔躍升邊界條件的烘房三維模型, 對青花椒堆積式烘干流場及均勻性具有較高的可靠性． 空載和加載時監測截面上基于模擬風速和實測風速的速度不均勻系數平均相對誤差分別為3.65%和2.09%．

2) 導流板參數對烘房內速度不均勻系數的影響從大到小依次為導流板角度、 導流板長度、 導流板位置． 以速度不均勻系數最小為目標, 導流板最優參數為導流板角度86.01°、 導流板位置924.75 mm、 導流板長度855.51 mm．

3) 導流板使烘房上部烘干區域形成2個子區間, 有利于烘干區域流場均勻性的提高． 與原始烘房相比, 優化烘房(空載)速度不均勻系數降低39.2%, 優化烘房(加載)速度不均勻系數降低20.1%, 烘房內流場均勻性得到較大提高．