徐靜安 吳文華 陳明良

技術講壇

第十五講基于UD+CFD耦合優(yōu)化的復混肥造粒塔進風口結構設計探討

徐靜安 吳文華 陳明良

0 引言

熔體高塔造粒過程主要是以熔融尿素為載體，與配比好的其他物料在熱熔狀態(tài)下攪拌混合，形成具有一定流動性的低溫共熔料漿，再通過旋轉差動噴頭將料漿以液流的形式噴出。液流在塔內自然下落過程中依靠表面張力成型，并與上升的空氣熱交換，冷卻凝固成球形顆粒。復混肥在造粒塔內的冷卻效果直接影響整個造粒過程，塔底復混肥黏壁、顆粒不圓整等均與塔內復混肥冷卻效果欠佳有關。而影響塔內復混肥冷卻過程的因素主要有進風口氣溫、復混肥顆粒粒徑、進出風口結構、塔高等。造粒塔的結構不僅影響著復混肥冷卻過程，也與造粒塔的造價成本息息相關。因此，在造粒塔設計階段，有必要對造粒塔結構（如進風口結構等）進行研究分析，嘗試找到合適的造粒塔結構。

迄今為止，國內已建成上百套造粒塔，但關于造粒塔結構對復混肥冷卻過程影響的研究甚少，尤其是對造粒塔透風窗結構方面的研究。國內不少學者對造粒塔內顆粒冷卻過程進行理論分析，將造粒塔內的空氣流動簡化成活塞流，在計算塔內風量時往往忽視了造粒塔進風口結構對塔內風量的影響，這勢必將降低計算結果的準確度。上海化工研究院長年致力于復混肥高塔造粒工藝研究工作，近年來采用計算流體動力學（CFD）軟件與現(xiàn)場測試數(shù)據(jù)相結合的方法建立了一套數(shù)學模型，該數(shù)學模型根據(jù)造粒塔的實際幾何結構進行建模，模型計算預測出的結果與現(xiàn)場測試數(shù)據(jù)的誤差在可接受范圍內。鑒于造粒塔設備龐大，塔高90～120 m，通過實驗及實測來掌握其性能工作量很大，本文在己建立的流體力學數(shù)學模型基礎上，通過耦合優(yōu)化方法改變幾何模型，初步探索研究造粒塔進風口結構對復混肥熔融造粒塔內風量的影響，以期為造粒塔幾何結構設計提供技術參考。

耦合優(yōu)化的邏輯框圖見本刊2017年第五期出版的第十三講“對基于‘試驗設計+模擬計算’耦合優(yōu)化研究、設計方法的推廣”。

1 初始模型計算方法

1.1 造粒塔初始結構

造粒塔內區(qū)域如圖1所示。造粒塔由混凝土筒體、進風口、出風口、旋轉差動噴頭和錐形接料斗等組成。造粒塔初始結構參數(shù)如表1所示。

圖1 造粒塔結構示意圖

1.2 網(wǎng)格劃分

造粒塔的模擬區(qū)域比較大，采用分區(qū)結構畫網(wǎng)格。將噴頭以上區(qū)域定義為氣流區(qū)，該區(qū)域網(wǎng)格為非結構化網(wǎng)格(Tet/Hyblrid)。將噴頭以下區(qū)域定義為顆粒區(qū)，該區(qū)域網(wǎng)格為結構化網(wǎng)格（Hex/Wedge）。另外，根據(jù)網(wǎng)格尺寸及數(shù)量對造粒塔進行網(wǎng)格獨立性檢驗，根據(jù)計算結果和計算時間選擇恰當?shù)木W(wǎng)格劃分方法。本文根據(jù)造粒塔尺寸共進行了三種水平的網(wǎng)格劃分，并對不同網(wǎng)格水平的的造粒塔模型分別進行模擬計算。鑒于網(wǎng)格數(shù)較大,常規(guī)計算機無法勝任，本次模擬計算工作均在上?；ぱ芯吭篒BM高性能計算機平臺上完成。具體網(wǎng)格數(shù)、模擬計算時間及模擬結果見表2。模擬結果表明，Level 2的網(wǎng)格數(shù)（即7543158）具有一定的合理性，在此基礎上增大網(wǎng)格數(shù)，模擬誤差減小不明顯，而計算時間增長，且殘差值增大甚至可能無法收斂?？紤]時間消耗，對造粒塔結構探索工作中選擇網(wǎng)格密度水平為Level 2。

表1 造粒塔初始結構參數(shù)

表2 網(wǎng)格獨立性檢驗結果

1.3 顆粒追蹤方式

本文數(shù)值模擬計算中，連續(xù)相均采用穩(wěn)態(tài)的計算方式，顆粒追蹤方式分別采用穩(wěn)態(tài)追蹤和非穩(wěn)態(tài)追蹤。顆粒穩(wěn)態(tài)追蹤時，“Number of Continuous PhaseIterations Per DPM Iteration”設置為20，即分散相每迭代一步連續(xù)相迭代20步。

顆粒非穩(wěn)態(tài)追蹤時，時間步長設置為0.1 s，時間步數(shù)為1，即每隔0.1 s顆粒噴射一次，每0.1 s對每個顆粒進行一輪包含一步的軌跡計算，每一輪顆粒軌跡計算中連續(xù)相迭代20步。

顆粒非穩(wěn)態(tài)追蹤與顆粒穩(wěn)態(tài)追蹤的模擬計算結果對比如表3所示。

表3 顆粒非穩(wěn)態(tài)追蹤模擬結果對比

顆粒非穩(wěn)態(tài)追蹤與穩(wěn)態(tài)追蹤的數(shù)值模擬結果非常接近，模擬結果與實測結果的相對誤差均在5%以內，其中顆粒非穩(wěn)態(tài)追蹤的模擬結果誤差相對較小。理論上，顆粒穩(wěn)態(tài)追蹤僅適合于顆粒稀少且顆粒對連續(xù)相影響較小的氣粒多相流的數(shù)值模擬，而顆粒非穩(wěn)態(tài)追蹤在模擬計算顆粒較密集、顆粒對連續(xù)相有一定影響的氣粒多相流領域更勝一籌。由于造粒高塔體積龐大,顆粒下降速度快、在塔內停留時間短，塔內氣粒多相流中顆粒所占體積分數(shù)較小，顆粒穩(wěn)態(tài)追蹤方式也適合于用來模擬計算塔內顆粒運動軌跡，最終體現(xiàn)在顆粒非穩(wěn)態(tài)追蹤與穩(wěn)態(tài)追蹤兩種方式模擬計算塔內氣粒多相流的模擬結果非常接近?？紤]到非穩(wěn)態(tài)計算耗時長且對計算機的性能要求高，本文在對造粒塔結構的模擬探索計算中均采用顆粒穩(wěn)態(tài)追蹤的方式。

1.4 物性參數(shù)及邊界條件

1.4.1 氣相定解條件

對于氣體入口與出口邊界條件，只需給定入口壓強p、溫度T、湍流強度I以及水力直徑Di。以其中一個模擬計算為例，氣相進口條件及參數(shù)值為：p= 101 325 Pa，速度方向為垂直于進風口截面，T=296.9 K，I=0.04，Di=3.4 m。氣相出口條件與參數(shù)值為：p=101 325 Pa，速度方向為垂直于出風口截面，氣流溫度T=296.9 K，I=0.05，Di=6.8 m。

1.4.2 顆粒定解條件

顆粒邊界條件主要包括顆粒進口條件和顆粒與壁面碰撞條件。

（1）顆粒進口條件

顆粒進口流量為10.56 kg/s，溫度為386.85 K (113.7℃)，假設固體顆粒除了比熱容外其他各種性質與料漿相同，進口物料熔液物性參數(shù)見表4。

顆粒從錐臺射出，速度大小為3.31 m/s，切向分量百分數(shù)0.999 989 1。錐臺注射的顆粒粒徑分布參考塔底顆粒粒徑分布（如表5所示）。

表4 復混肥物性參數(shù)

表5 塔底顆粒粒徑分布

模擬計算過程中將顆粒粒徑分布近似成R-R分布，遵從R-R分布的顆粒群中，顆粒直徑d與尺寸小于d的顆粒所占的質量分數(shù)Rd有如下關系：

式中，d為特征尺寸，其含義為直徑小于d的顆粒質量分數(shù)Rd=1-e-1=0.632 1；n為尺寸分布指數(shù)，該參數(shù)值越大，顆粒尺寸分布越均勻。數(shù)值模擬過程中，顆粒群的特征尺寸d為2.48 mm，分布指數(shù)n為3.5。

復混肥顆粒在塔內的降落過程還伴隨著凝固放熱過程。顆粒從外往里開始凝固，因此存在開始凝固溫度Tf和完全凝固時的平均溫度Ts。近似地認為顆粒溫度高于Tf時為液態(tài)，溫度低于Ts時為固態(tài)，凝固潛熱采用修改比熱容的方式加入。具體修改方式如下：

式中，C1為液態(tài)料漿熱容，CS為固態(tài)顆粒熱容，L為凝固潛熱。

復混肥配方為20-5-20（硫基），該種復混肥的開始凝固溫度Tf=86.1℃，完全凝固時的平均溫度Ts=83.7℃，凝固潛熱L=87.34 kJ/kg。

（2）顆粒出口條件

簡化模型中將顆粒出口面設置為“壁面”，顆粒與出口面的碰撞條件取“逃逸”，也就是顆粒到達出口面時計算終止。

（3）顆粒碰撞條件

液滴進入造粒塔內，與空氣逆流接觸實現(xiàn)熱量交換，由于現(xiàn)有工藝條件下料漿液滴較少黏壁，固體肥料碰壁后假設成彈性碰撞，因此顆粒與壁面的碰撞條件取“彈性碰撞”。

1.5 計算原理

采用分散相模型進行計算，控制方程包括連續(xù)相方程、分散相方程和氣粒兩相間熱量傳遞方程。連續(xù)相求解過程為穩(wěn)態(tài)迭代計算，顆粒追蹤方式采用穩(wěn)態(tài)追蹤。所謂穩(wěn)態(tài)方式是指每隔若干個連續(xù)相流場迭代步，在當前流場狀態(tài)下，逐個對每個顆粒進行從起始位置直到其終了的軌跡積分計算及源項計算。穩(wěn)態(tài)方式得到某一時刻連續(xù)相流場條件下，在一系列積分時間步的顆粒狀態(tài)，一系列顆粒位置可連成運動軌跡線。

熱壓自然通風主要是依靠浮升力作用的一種通風方式，自然通風的數(shù)值模擬研究方法一般為擴大計算域，在建筑的外圍建立一較大的模型形成無限遠邊界，在外圍入口處設置為速度入口邊界條件來考慮側風對其的影響。不考慮側風的影響時，可以將速度邊界的速度值設得很小，或者改為壓力邊界條件，或者直接將建筑的進出口處均設為壓力入口邊界條件。綜合考慮后，決定將造粒塔進出口邊界條件均設為壓力邊界條件，內部空氣密度也采用波希涅斯克(Boussinesq)假設。

離散方式采用有限差分法，連續(xù)相方程采用SIMPLE算法迭代求解，分散相建立在連續(xù)相初始求解值上，并與連續(xù)相進行耦合求解。

2 造粒塔初始條件模擬結果分析

本文針對兩座工程塔的三次不同工藝參數(shù)下的造粒過程分別進行模擬計算，具體模擬結果及工程標定誤差如表6所示。其中“模擬1”與“模擬2”針對的造粒塔進風口結構：n=9，s=74.46 m2、θ=139.81°；“模擬3”針對的造粒塔進風口結構：n=14，s=90.72 m2，θ=160.51°。表6顯示，模擬結果與工程化高塔實測結果的相對誤差均在10%以內，說明建立的數(shù)值模型可靠，得出的數(shù)值模擬結果具有一定的可信度。

表6 數(shù)值模擬結果對比

具體針對“模擬1”文件進行造粒塔內氣速氣溫云圖分析，如圖2所示。并對造粒塔內顆粒分布、運動降溫規(guī)律等進行分析，如圖5與圖6所示。造粒塔外面的冷空氣在壓差的作用下，從9個進風口進入造粒塔，9個不同氣速的氣流在塔中心匯流形成渦流。從圖2中可以明顯看出塔中心的渦流，使得粒子發(fā)生錯流，軌跡發(fā)生變化，部分粒子可能從氣速小的窗口逃逸，現(xiàn)場測試中也發(fā)現(xiàn)通風窗口落有不少粒子。氣體匯流后沿著塔軸向上升，上升過程中與下降的復混肥顆粒進行熱量交換與動量交換。

圖2 塔進風口處橫截面上氣速矢量圖

圖3和圖4顯示的是造粒塔內xz截面溫度云圖、速度云圖、溫度等值線和速度等值線。從圖3中可以看出溫度等值曲線呈現(xiàn)規(guī)則的W形，出現(xiàn)兩處低谷一處高峰。這與造粒塔顆粒分布規(guī)律有關，如圖5所示，顆粒主要集中在塔半徑4～7 m區(qū)域內，顆粒分布規(guī)律是由噴頭噴出顆粒的粒徑分布規(guī)律導致的，粒子的粒徑主要集中在2～3.36 mm。塔半徑4～7 m區(qū)域內粒子比較多，導致該區(qū)域空氣溫度上升得比較快，升高到相同溫度所需上升的高度比起他區(qū)域要低，最終導致軸截面的溫度等值線中出現(xiàn)兩處低谷。塔中心粒子很少，熱交換不明顯，氣溫升高到相同溫度所需升高的高度比其他區(qū)域要高，最終導致軸截面的溫度等值線中出現(xiàn)高峰。氣體上升過程中，溫度也在上升，在噴頭區(qū)域溫度梯度明顯變大。這是因為噴頭區(qū)域下降的是高溫熔融復混肥并且顆粒比較密集，復混肥與空氣的溫差比較大，傳熱推動力大。噴頭區(qū)域的熱空氣繼續(xù)向上運動，由于漸漸遠離熱源以及與塔壁附近的冷空氣換熱，氣溫逐漸下降。

圖3 塔xz軸截面上氣溫云圖

圖5 距進風口35 m處橫截面上顆粒分布圖

從圖4中可以看出進風口氣速最大，空氣進入造粒塔后由于截面迅速增大導致速度迅速減小，空氣上升過程中沿軸向速度變化不是很大。進入噴頭區(qū)域由于熱壓增大，因此氣體明顯提速，進入熱氣通道后橫截面減小，氣速再次提高。塔半徑4～7 m處氣速明顯較小，這與粒子集中在該區(qū)域，空氣受到較大的阻力有關。

圖6 塔內顆粒降溫規(guī)律圖

顆粒在造粒塔內下降并與上升空氣發(fā)生對流換熱，顆粒溫度逐漸降低，依次經(jīng)歷液態(tài)冷卻過程、凝固過程、固態(tài)冷卻過程。從圖6可以看出2.0 mm粒徑顆粒下降5 m左右進入凝固過程，下降35 m左右進入固態(tài)冷卻過程；2.5 mm粒徑顆粒下降10 m左右進入凝固過程，下降50 m左右進入固態(tài)冷卻過程；而3.9 mm與4.1 mm粒徑顆粒下降25 m左右進入凝固狀態(tài)，落至塔底時依舊處于凝固狀態(tài)，顆粒溫度最終為86.2℃。不同粒徑顆粒的傳熱系數(shù)與停留時間的差異最終導致了顆粒降溫曲線之間的差異。大粒徑顆粒對流傳熱系數(shù)與停留時間均小于小粒徑顆粒的對流傳熱系數(shù)與停留時間，降落至塔底時凝固潛熱尚未完全釋放出來。

3 進風口結構的耦合優(yōu)化探索

3.1 均勻試驗設計

本文主要是在已有的被校核的模擬計算基礎上，通過改變進風口結構分析研究進風口結構對自然通風造粒塔內風量V的影響，其中造粒塔結構參數(shù)及工藝參數(shù)如表7所示?？紤]到影響塔內風量的因素較多，只選取對塔內風量影響較大的3個因素。即通風窗個數(shù)n、進風口總面積s以及通風窗對應的總圓心角θ。各因素選取水平參照均勻設計表U7*(74)，各因素具體取值如表8所示。

表7 結構參數(shù)及工藝參數(shù)

表8 因素水平表

其中造粒塔具體進風口結構如圖7所示，造粒塔底部徑向分布著若干個矩形通風窗，通風窗個數(shù)n一般為6～18個，進風口總面積s一般為66～96m2，通風窗對應的總圓心角θ一般為100～180°。

圖7 某造粒塔進風口結構布置圖

3.2 均勻設計結果分析

對于高約l00 m的工程造粒塔，現(xiàn)場實測標定一次工作量相當大，涉及到結構尺寸的變動改造幾乎不可能。而基于耦合優(yōu)化方法的均勻設計方案，開展對應的數(shù)值模擬計算，以造粒塔內風量V為分析指標可以得到相應的技術信息，具體結果見表9。

表9 模擬試驗結果

通過數(shù)據(jù)處理系統(tǒng)（DPS）軟件對模擬試驗結果進行二次多項式逐步回歸統(tǒng)計建模，計算結果如下所示：

計算結果表明各項檢驗顯著，模型擬合效果在可接受范圍之內。t檢驗顯示在實驗范圍內對V的影響次序：s×θ＞n×θ。數(shù)學統(tǒng)計模型顯示造粒塔內風量隨著通風窗個數(shù)增大而減小，隨著通風窗總面積與總圓心角的增大而增大。此外,數(shù)學模型預測的最優(yōu)化進風口結構為：n=6，s=96 m2，θ=180°。

4 耦合優(yōu)化結果驗證

統(tǒng)計模型耦合優(yōu)化預測的進風口結構為n=6，s=96 m2，θ=180°，對造粒塔重新進行幾何建模，并結合對應的工藝參數(shù)進行模擬計算，模擬結果顯示優(yōu)化結構下造粒塔風量為152.35 m3/s，模擬計算值與統(tǒng)計數(shù)學模型預測值的相對誤差為1.05%。直接根據(jù)優(yōu)化進風口結構建造造粒塔，造價上千萬。鑒于直接工程驗證優(yōu)化結果成本太大，本文選擇對已有的山東某造粒塔進行實測標定，驗證耦合優(yōu)化模型預報的質量。山東造粒塔進風口結構為：n=9，s=74.46 m2，θ=139.81°，生產(chǎn)工藝參數(shù)如表7中所示。根據(jù)數(shù)學模型V=141.649 736 6-0.003 530 702 055×n×θ+ 0.000 932 333 124 6×s×θ計算，計算結果如表10所示，結果顯示模型預測值與實測標定值的相對誤差為0.975 4%，誤差在可接受范圍內，說明回歸的耦合優(yōu)化模型具有顯著的統(tǒng)計意義。

5 結論

（1）介紹了復混肥造粒塔CFD數(shù)值模擬計算原理與過程，并分析研究了造粒塔內空氣溫度、速度云圖以及復混肥顆粒在塔內的運動規(guī)律。分析結果顯示造粒塔內軸截面氣溫等值線呈現(xiàn)W形狀，氣速等值線呈現(xiàn)M形狀，在進風口處造粒塔中心與通風窗的上端均出現(xiàn)明顯的渦流。顆粒在造粒塔內降落過程中主要集中在塔半徑4～7 m的范圍內，依次經(jīng)歷液態(tài)冷卻過程、凝固過程、固態(tài)冷卻過程，顆粒溫度隨路程的變化曲線出現(xiàn)兩個明顯的轉拆點。

表10 模型預報驗證結果

（2）基于均勻設計法開展數(shù)值模擬的耦合優(yōu)化研究，探索造粒塔進風口結構對塔內風量的影響。通過對模擬結果進行二次多項式回歸分析，得出數(shù)學統(tǒng)計模型，并最終確定了試驗范圍內的最優(yōu)化進風口結構：n=6，s=96 m2，θ=180°。

（3）結合工程實測標定，對耦合優(yōu)化模型預報質量進行驗證，結果顯示模型預測值與實測標定值的相對誤差在可接受范圍內，回歸的模型具有顯著的統(tǒng)計意義，可作為工程設計技術基礎。

后記：

最近，閱讀了《中國人才》2015年9月刊登的一篇文章——各國加強科技人才隊伍建設特點與趨勢（上）。各國均在創(chuàng)新人才培養(yǎng)方面強化科學、技術、工程和數(shù)學(STEM)教育，借此推動國家創(chuàng)新和經(jīng)濟的發(fā)展。英國政府在2011年—2015年實施支持STEM教育計劃，美國2013年公布了聯(lián)邦STEM教育五年戰(zhàn)略計劃，加拿大政府2014年12月啟動STEM的人才培養(yǎng)計劃，韓國在2011年—2015年強化科學、技術、工程、藝術、數(shù)學(STEAM)人才的系統(tǒng)培養(yǎng)。

多年來，上?；ぱ芯吭褐С趾屯苿訑?shù)學工具、數(shù)字化技術在應用基礎研究、技術開發(fā)、工程化轉化工作中的學習、探索、應用。2010級碩士吳文華在讀期間及留院工作后，結合高塔造粒項目，不斷探索UD+CFD耦合優(yōu)化技術，并在院部配置的IBM高性能計算機上重新做了計算，不斷完善。

“沉舟側畔千帆過，病樹前頭萬木春”，青年同仁們，努力??！

徐靜安男教授原上?；ぱ芯吭涸洪L長期從事化工機械、化學工程、化工工藝開發(fā)、化工試驗設計與數(shù)據(jù)處理方面的研究