氣固流化床啟動階段擋板內構件受力特性的CFD-DEM 模擬

李鐵男，趙碧丹，趙鵬，張永民，王軍武，4

（1 中國石油大學（北京）重質油國家重點實驗室，北京 102249； 2 中國科學院過程工程研究所多相復雜系統國家重點實驗室，北京 100190； 3 中國科學院大學化學工程學院，北京 100049；4中國科學院綠色過程制造創新研究院，北京 100190）

引 言

流化床反應器具有良好的傳質傳熱性能，且能夠較好輸送固體顆粒物料，因此廣泛用于工業生產過程[1]。反應器內設置合適的內構件可強化傳質與傳熱并提高反應器性能[2-3]?；谝酝芯拷涷?，采用傾斜葉片的流化床內構件破碎氣泡和改善氣固接觸的效果更好，而且還可以通過調整結構參數有效調節氣固相的停留時間分布，具有比垂直構件和具有豎直葉片的網狀格柵更好的反應強化效果，因此此類內構件在工業中應用最廣泛[3]，典型的例子如百葉窗格柵、脊型和塔型內構件等。由于工業流化床內存在復雜的動態兩相流動以及不同操作狀態之間的切換，內構件會受到不同形式和強度的作用力，如果內構件結構設計不合理，很可能會致使內構件發生損傷甚至被破壞，進而導致工業裝置出現故障和經濟損失。因此，為了保障流化床內構件的長周期可靠性，系統掌握內構件在流化床內的受力特性尤為重要[4]。

早期針對流化床中內構件受力特性的研究主要采用實驗手段，大多數關注的是正常流化狀態下B 類顆粒流化床內圓形水平換熱管的受力特性[5-9]，其中Grace等[7-8]的研究為揭示正常流化狀態下圓管受力的機理奠定了良好基礎。近年來，Zhang等[4,10-13]系統研究了化工領域A 類顆粒流化床反應器中經常使用的斜片擋板內構件的受力特性，發現啟動階段流化床層內部構件會受到一個很大的向上載荷脈沖，其峰值可達正常流化狀態下內構件受到的平均載荷值的數倍甚至一個數量級以上[11]，這對工業裝置內構件可靠性是一個巨大的潛在威脅，很可能是導致工業裝置內構件損壞的一種新機理。

對于冷態實驗室條件下一些小型簡單結構內構件，實驗方法可以獲得較為準確的受力特性數據，但對于工業裝置高溫高壓的苛刻條件下的大型復雜結構內構件，其在流化床中的受力特性預測則難以使用實驗方法，因此，使用近年來日益成熟的計算流體力學模擬方法研究流化床中內構件的受力是一種有益的嘗試。Higashida等[14]采用虛擬顆粒方法（FPM 方法）模擬研究了三維鼓泡床中懸浮球體受到的動態垂直力，并將模擬計算結果與帶有拉格朗日傳感器系統的懸浮球體所測得的實驗值相比較。研究發現，球體在裝置內床層流化過程中，其動態垂直力的波動與氣泡的運動密切相關，球體受到的合力由流體作用力與顆粒碰撞作用力兩部分組成。Yan 等[15]采用雙流體模型模擬了氣固兩相體系內不同形狀的大顆粒受力演化規律。模擬結果顯示：隨著表觀氣速的增加，大顆粒受到的作用力也相應增大。Nagahashi 等[16]采用FPM 模擬方法，對二維流化床中單個氣泡經過水平圓管的受力機理進行了模擬研究，并將模擬結果與實驗結果相比較。結果顯示：模擬與實驗中氣泡的運動軌跡吻合較好，實驗拍攝氣泡經過圓管過程中，氣泡尾流攜帶的顆粒撞擊圓管會使圓管受到一個脈沖峰值。

本研究以啟動階段流化床層內部構件的非正常受力現象為研究對象，首先基于CFD-DEM 數值模擬方法，建立流化床擋板內構件表面載荷強度的統計方法，并通過與前期流化床啟動階段斜片擋板內構件受力特性的實驗數據[11]進行對照，以驗證該載荷統計計算方法的正確性及合理性。進一步，以該方法研究不同參數（碰撞模型參數、表觀氣速、顆粒粒徑）對啟動階段斜片擋板內構件受力特性的影響規律，旨在為下一步揭示流化床啟動階段內構件非正常受力現象的內在機理奠定良好基礎，也為開發工業流化床內構件受力特性的預測工具奠定理論基礎。

1 擋板受力統計計算方法

氣固擋板流化床啟動階段的本質是密相顆粒料層由未流化的固定床狀態逐漸轉變為流化床狀態的一個過渡過程，相比其他多相流模擬方法，CFD-DEM 方法更適宜于準確描述這一過渡過程。因此，本研究將基于CFD-DEM 方法，建立內構件表面載荷分布的統計計算方法，分析內構件表面受到的動態載荷信號。從統計力學角度分析，氣固流化床密相床層中的內構件受力特性是氣相分子與離散顆粒對內構件壁面作用特性的動態表征。CFDDEM 模擬方法可以獲得顆粒速度、位置和相互作用力(顆粒-顆粒、顆粒-壁面)等微觀顆粒相尺度信息，但是無法直接獲得工程上關注的內構件表面載荷強度及分布等宏觀物理量信息。

氣體壓力的物理本質為氣固兩相流體系內構件受力求解提供一定參考。內構件氣固兩相流體系內不同尺度下，從統計物理學可知固體顆粒與內構件表面作用的動態變化過程和氣體分子與壁面作用的動態過程相類似。此外，宏觀連續場信息通常都是借助顆粒流系統中每個離散顆粒和分子動力學中單個原子的作用力、位置和速度等微觀信息而獲得，故內構件表面受到的顆粒碰撞作用力特性可參考氣體對壁面作用形成的壓力過程及其動態作用力特性。在標準分子動力學模擬中，基于維里定理求解計算表面壓力。對于非均勻系統，壓力應為張量形式，且是與位置相關的函數[17]。其中，壓力張量的公式包含兩部分：第一部分稱為動壓力，代表顆粒的運動對壓力的貢獻量；第二部分稱為位形壓力，代表粒子間的相互作用。

氣固兩相流體系內構件受力載荷的求解計算需要了解顆粒流體系內應力張量的組成形式。連續介質場通常需要由離散粒子數據構造，這些離散數據是每個原子或顆粒的位置、速度和相互作用力。在顆粒系統中，顆粒動理學理論是在Chapman等[18]的氣體動理學理論的基礎上發展起來的。一般而言，封閉固相應力的方法主要有兩類：一類是顆粒動理學理論(KTGF)；另一類是基于CFD-DEM 方法的數值統計分析。KTGF 的推導方法主要分為以下幾種。(1) 考慮非彈性碰撞及顆粒體積影響。Gidaspow[19]基于Chapman等[18]建立的氣體動理論，進而考慮非彈性碰撞影響和顆粒體積影響等不可忽略因素，建立起經典的顆粒動理論。(2)考慮固相處于近平衡穩態情況。Rao 等[20]進一步考慮固相處于近平衡穩態時，通過完善顆粒速度分布函數，導出更合理的固相本構關系。(3)考慮體系內介尺度結構廣泛存在情況。Zhao等[21]考慮氣固兩相流系統中團聚物等介尺度結構廣泛存在，建立了適用性更廣且更準確的多尺度動理論。而對于另一類基于CFD-DEM 方法的數值計算，可通過捕捉顆粒系統中單顆粒的運動特性，選擇合理的統計平均方式，獲得合理的固相本構關系。顆粒系統中應力張量的求解方法主要分為三種。（1）體積平均法。Babic[22]最早提出了體積平均法，而后Zhu 等[23]在此基礎上進一步發展，將顆粒質量、速度和相互作用力等微觀性質與體系密度、應力和偶應力等宏觀特性聯系起來。（2）平面法。Mehrabadi 等[24]采用虛功原理應用平面法推導出應力張量表達式，并證明平面法推導出的應力張量形式與體積平均法等效。（3）粗?；椒╗25]?；诩僭O顆粒間為二元碰撞，每一對顆粒都只有單一的一個接觸點，碰撞接觸區域被一個接觸點取代，且碰撞不是瞬時的，將離散數據體積平均化并用連續性描述。此方法推導應力張量表達式不需要假設顆粒是剛性的或球形的。綜上所述，顆粒體系內的固相應力張量的組成形式主要包含兩個部分：(1) 動應力(kinetic stress)，由顆粒體系內局部體積中顆粒質心的速度相對于局部體速度的脈動產生；(2) 接觸應力(contact stress)，由體系內顆粒之間的直接接觸產生。

當顆粒體系內有擋板內構件存在時，內構件邊壁與顆粒間的碰撞也會對應力張量有貢獻作用。為了建立本研究中擋板內構件表面載荷的統計計算公式，需要了解內構件的存在對應力張量的貢獻形式。Weinhart等[26]采用Goldhirsch[25]的體積平均方法導出了邊界存在所貢獻的邊界應力形式。研究發現，當體系內有壁面邊界存在下，應力張量除了由顆粒速度波動產生的動應力和流動顆粒間碰撞產生的碰撞應力外，還有邊壁與顆粒碰撞的碰撞貢獻。根據牛頓第三定律可知，壁面邊界對固相顆粒的碰撞作用與顆粒施加于壁面的碰撞作用數值上相等。因此，統計擋板表面受到顆粒的碰撞載荷形式與邊界碰撞應力的統計形式類似。

統計計算密相床層中的擋板內構件表面受到顆粒的碰撞力，旨在將微觀的顆粒壁面碰撞信息與宏觀的擋板表面受力信息聯系起來，并將信息在時間與空間上統計平均得到內構件表面受到顆粒碰撞應力及其分布特性。本研究的擋板內構件形式為斜片擋板，其受到顆粒的碰撞作用主要表現在與氣體分布板通入氣體的流向相垂直的上、下兩表面所受的顆粒碰撞力。

根據柯西基本定律，擋板內構件表面上單位面積受到的顆粒碰撞力與擋板表面面積平均的碰撞應力張量關系如式（1）所示。

式中，v為單位向量；σ為碰撞應力張量；T為碰撞力。

基于擋板內構件表面受到的碰撞載荷與固相應力中內構件邊壁的貢獻數值上相等，故需要找到一個合理的平均化方法將顆粒的微觀碰撞信息轉化為擋板表面受到的宏觀應力信息。對體系內邊壁存在的影響處理，Ries 等[27]提出了鏡像系統(the mirrored system)的方法。此方法的思想是通過類似將體系內邊壁一側的顆粒關于邊壁鏡像對稱，來代替邊壁存在對固相應力的影響。因此，基于已知擋板內構件的表面面積，為確定擋板表面面積平均的碰撞應力張量σ，可參考此方法處理內構件表面。

由于本研究所采用的顆粒物料為單一粒徑的顆粒(帶粒徑分布的也可以以同樣的方法統計分析，但是此時兩顆粒質心間位移連線的長度l不是常數)，為了更合理有效地統計計算擋板表面受到的碰撞應力，將擋板表面緊貼著的一層顆粒關于擋板表面鏡像對稱，此時與顆粒相互作用的擋板表面等效為一層粒徑與床內顆粒相同且平鋪面積與擋板表面面積相同的平面顆粒。因此，擋板表面受到顆粒的碰撞作用等效為單層平面顆粒與虛擬顆粒之間的碰撞作用，可認為單層平面顆粒與碰撞顆粒組成的顆粒群間存在一個法向量與擋板表面法向量相一致的“虛擬平面”，如圖1所示。

圖1 密相床層中擋板表面與顆粒接觸示意圖(a)；等效為“虛擬平面”后顆粒碰撞接觸圖示(b)Fig.1 Schematic diagram of contact between the baffle and particles in dense bed(a)and contact between an imaginary plane and particles in a dense bed(b)

水平放置的擋板主要以上、下兩個表面為碰撞受力面，故在此以擋板上表面為例（下表面處理方法一致）。單位面積的“虛擬平面”上存在足夠多數量的顆粒間碰撞，定義T為單位“虛擬平面”面積上碰撞顆粒對單層平面顆粒施加的碰撞力，則T如式（2）所示。

式中，Fi為單位面積“虛擬平面”的上、下第i對顆粒的碰撞接觸力；Nc為單位面積上的碰撞接觸點數。

Mehrabadi 等[24,28]通過虛功原理證明了單位“虛擬平面”面積上顆粒碰撞力的統計求和可寫成統計平均的顆粒簇體積V上的求和，此采樣顆粒簇由與擋板上表面接觸的真實顆粒和鏡像虛擬擋板顆粒組成?！疤摂M平面”上面積平均應力張量σ與采樣體積內顆粒間相互作用產生的碰撞應力張量σ′相等[29-30]，即

因此，由邊壁表面受到的顆粒碰撞作用貢獻形式可知，統計平均的顆粒簇體積中顆粒間碰撞應力張量σ的統計計算公式為

式中，上角標w 代表壁面;下角標α和β表示方向分量。

本研究中，擋板上、下表面的表面積均為S；顆粒粒徑為dp；一對顆粒的每次碰撞接觸中，兩顆粒質心間位移連線的長度用l表示，即l的大小為顆粒質心到碰撞接觸點間位移長度的2 倍；F為兩顆粒間的碰撞接觸力；Vw為擋板表面相對應的均勻統計顆粒簇體積，如圖1(b)中虛線框所示，其值大小為以擋板表面積為底，顆粒粒徑為高的立方體體積，即Vw=Sdp；N為某一時刻擋板表面的總碰撞接觸點數。因此，擋板表面受到顆粒間碰撞應力張量σ統計計算公式為

從式（5）可知，在數值模擬計算程序中，每一時間步求得擋板表面受到顆粒碰撞的應力張量表達式為

式中，每個分量表示擋板表面受到的均布載荷；σ33為垂直擋板表面方向的碰撞均布載荷，為主要關注參數。運用此公式統計分析CFD-DEM 的模擬結果，即可分析出顆粒與壁面作用的應力分布特性。

2 模擬的幾何模型及物理參數設置

本研究選取Liu 等[11]的三維方形冷模流化床實驗裝置（圖2）作為模擬研究對象。裝置的主體由方形流化床床體、預分布器、氣體分布板、單個斜片擋板和旋風分離器等組成。方形流化床的截面面積為0.3 m×0.3 m，總高度為5 m。單個斜片由單個板條和兩個固定底座組成，固定底座通過螺栓將板條的兩端固定在流化床的床壁上，使得板條可穩定置于密相床層中。其中，擋板板條的長為0.3 m、寬為0.05 m、厚度為0.008 m，且板條的長度可剛好橫跨床層截面，與流化床的截面邊長一致。單個斜片擋板水平安置在流化床邊壁的中心位置處，安裝高度距流化床底部分布板0.5 m。

圖2 三維方形冷模流化床實驗裝置示意圖[11]Fig.2 Schematic diagram of the cold three-dimensional fluidized bed with a square cross-section[11]

為了便于計算，本研究將上述實驗裝置進行簡化處理，省略了流化床頂部的旋風分離器和底部的氣體分配室。方形流化床幾何模型的長度和寬度與實驗等同，即長和寬均為0.3 m。實際實驗中流化床內床層膨脹高度最高膨脹至2 m 左右，考慮到節省計算資源和提高計算效率，將流化床的床層稀相空間的高度縮減，故方形流化床的模型總高度定為2.2 m。模型中，擋板內構件的安裝位置與實驗相一致，安裝在距床層底部分布板0.5 m 的流化床壁面中心位置處，且擋板長度保持不變，仍為0.3 m。為簡化網格劃分并成功生成均勻六面體網格(六面體網格尺寸為0.01 m×0.01 m×0.01 m)，故將單個斜片擋板內構件的寬度調整為0.06 m，厚度定為0.01 m。而對于本研究中的無擋板自由床而言，除略去擋板流化床內斜片擋板外，其余幾何結構均與擋板流化床一致。最終得到模擬采用的方形流化床幾何結構如圖3所示。

圖3 模擬采用的流化床幾何結構Fig.3 The geometry of fluidized bed for simulation

本文選用劉對平[4]實驗中的粒徑為595 μm的B類非球形石英砂顆粒作為模擬體系內的固體介質，密度等相關物料基本性質與實驗測量值一致。實驗前期未曾測量石英砂顆粒的球形度信息，故根據測量顆粒的起始流化氣速等參數信息，采用Hua等[31]提出的方法求解估算顆粒的球形度，最終確定模擬采用的顆粒球形度為0.86。流化床中床層物料的初始自由堆積高度為1 m。模擬體系內流體采用常溫常壓空氣（25℃、1.2 kg/m3、1.8×10-5Pa·s）。

3 數學模型及模擬設置

3.1 數學模型

本研究所采用的粗?；疌FD-DEM 方法是基于Lu 等[32]提出的EMMS-DPM (energy-minimization multi-scale-discrete particle method，基于能量最小多尺度離散顆粒)方法。故本研究中的粗顆粒質量mCGP=k3mp，其中，k為粗?；剩╧=dCGP/dp），mp為真實顆粒的質量。顆粒-顆粒碰撞作用與顆粒-壁面碰撞作用的計算使用Peng 等[33]提出的方法，采用不考慮歷史碰撞影響的Hooke（胡克線彈性碰撞）模型求解計算。計算顆粒-顆粒碰撞與顆粒-壁面碰撞之間的相互作用力時，盡管非球形顆粒與顆粒和壁面之間碰撞接觸力的計算，可以通過組合顆粒球元構建非球形顆粒幾何模型的方法實現[34-36]，但將大幅增加計算量。因此，本文在計算顆粒-顆粒以及顆粒-壁面間相互作用力時仍然將顆粒視為球形顆粒。顯然，這是一種提高計算精度與保證合理計算量間采取的折中方法，研究表明這種方法可以獲得較為合理的模擬結果[37-39]。氣固相間曳力的計算采用Hua 等[31]提出的非球形顆粒的曳力模型求解。在Ergun 曳力關聯式中引入球形度，并利用Ganser非球形顆粒曳力模型[40]來計算Wen 和Yu 曳力關聯式中的單顆粒曳力系數Cd。

氣固系統DEM 的控制方程組分別為氣相、固相的質量守恒和動量守恒方程：

式中，ρg、ug和p分別為流體密度、速度和壓強；t為時間；g為重力加速度；εg為網格空隙率，即計算網格內流體的體積分數；NCGP,cell為計算網格中粗顆粒數目；Vcell為當前計算流體網格的體積；Fdrag,i為計算網格內顆粒i與氣流之間的曳力交換；τg為牛頓流體黏度應力張量；μg為流體的剪切黏度；I為單位張量；uCGP,i為粗顆粒平動速度；dCGP,i為粗顆粒粒徑；N為某一時刻與顆粒i相互作用的顆?？倲?ICGP,i為粗顆粒的轉動慣量；ωCGP,i為粗顆粒角速度；分別為法向和切向碰撞接觸力；R為顆粒質心到碰撞接觸點的位移矢量；μt為滾動摩擦因數；RCGP,i為粗顆粒的半徑。

鑒于氣固流化床內顆粒之間的碰撞接觸較為頻繁，本研究采用彈簧-阻尼模型對粗顆粒之間的碰撞過程進行求解計算[33]，顆粒所受的碰撞接觸力等于法向接觸力與切向接觸力之和。

本研究中固相顆粒為非球形的石英砂顆粒，故采用Hua 等[31]提出的非球形顆粒的曳力模型來確定本研究系統內氣固相間曳力。其中具體的曳力系數可表示為[41]

式中，εp為網格中顆粒的體積分數；Ψ為顆粒的球形度；mi為粗顆粒質量；ui為粗顆粒速度；k為粗?；?；Rep為顆粒Reynolds 數；K1和K2分別為顆粒斯托克斯形狀因子與牛頓形狀因子；dA為顆粒等投影面積球當量直徑；dV為顆粒等體積球當量直徑，且dA和dV兩參數難以通過實驗測得，故假設dA=dV=dsauter為顆粒Sauter 平均直徑；D為流化床的水力直徑；εi為網格內i類顆粒的體積分數；di為顆粒的粒徑。由于本研究體系內顆粒為單一粒徑顆粒，故di=dave，εi=εp。

3.2 模擬設置

本研究所需方形氣固擋板流化床與無擋板自由床的三維幾何結構由ANSYS ICEM 軟件進行建立，模擬參數如表1 所示。采用均勻的六面體結構化網格進行劃分繪制完成對整個計算床體的網格表征，而后將網格導入開源計算流體力學軟件OpenFOAM 進行求解計算。流體計算網格用于求解CFD 中計算流體的空隙率及其他物理量，對保證數值計算的精確性和穩定性具有重要意義。流體計算網格的尺寸需足夠精細方可獲得充足的尺度信息[41-42]，且為保證計算資源與時間成本的合理性，故取邊長為10 mm的六面體網格作為后續模擬研究的基礎，網格尺寸約為真實顆粒粒徑的15 倍，亦是粗?；蟠诸w粒的2～5 倍，故可保證用顆粒中心算法可合理統計計算網格中的空隙率及固相速度等物理量。分別在幾何模型的X、Y和Z三個方向上繪制30×30×220 個均勻流體計算網格。此外，本文所采用的固相顆粒物料石英砂顆粒為Geldart B 類顆粒，劃分的流體計算網格足以獲得網格無關性的模擬結果[43]。

表1 擋板床與自由床數值模擬參數Table 1 Simulation parameters for the baffled and free fluidized bed

圖3 顯示了本文數值模擬的邊界條件設置情況。流化床床層底部作為氣相流體的入口，將其設置為速度進口邊界條件，底部采用均勻進氣，且速度大小依據操作氣速設定。氣相出口位于反應器幾何模型的頂部，并設置為壓力出口邊界條件，且壓力值與大氣壓相同。為模擬前期實驗中二級旋風分離器及濾袋的作用，在反應器幾何模型頂部設定不允許顆粒逸出，以保證模擬體系內的顆粒存量恒定。擋板內構件設置為與流化床壁面一致，設定為無滑移的邊壁。為確保計算資源與時間成本的合理性，本研究選取模擬離散顆粒的粗?；蔾=5，且流體計算網格的尺寸是粗?；蟠诸w粒的2～5倍，可保證模擬計算的準確性[42,44]。根據流化床內真實顆粒數目確定的粗?；w粒數為3459863個。

4 模擬結果與討論

4.1 碰撞受力統計計算公式驗證

擋板表面法向應力分量σ33在顆粒與壁面的碰撞作用中占絕對主導地位，其統計結果的數量級遠遠大于其余8 個分量。因此，擋板表面法向應力分量σ33可作為評價擋板表面受顆粒碰撞作用大小的參數，其余8 個分量可忽略。此外，由前面分析可知，擋板表面法向應力分量σ33為擋板上、下表面受到顆粒碰撞的均布載荷差，即擋板受到的顆粒碰撞載荷為擋板下表面受到的碰撞均布載荷與上表面受到的碰撞載荷的差值，即

式中，q為擋板受到的顆粒碰撞均布載荷，方向豎直向上為擋板下表面受到的顆粒碰撞均布載荷；為擋板上表面受到的顆粒碰撞均布載荷。

從圖4 中可以看出，流化床啟動階段內構件所受碰撞載荷的CFD-DEM 模擬結果與實驗結果隨時間演變的趨勢大體上一致，數值模擬統計的擋板受均布載荷變化曲線顯示在流化氣體剛通入的啟動階段會產生一個較大的脈沖峰值，而后擋板內構件所受的載荷減小并進一步波動變化。與實驗值相比，啟動階段的第一個波峰及波谷值出現的時間更早，這說明應用本文所建立的載荷統計計算公式可以定性地成功模擬統計出床層由固定床向流化床轉變過程中內構件受到非正常受力特性這一現象。此外，圖中CFD-DEM 數值模擬結果的脈沖峰值與實驗值大小基本一致，說明模擬可以(半)定量復現實驗中觀測到的氣固流化床反應器中密相床層流態轉變時的受力特性。最后，在后續穩定階段載荷波峰的出現與氣泡尾流內顆粒與壁面作用緊密相關。綜上所述，擋板內構件受顆粒碰撞載荷統計計算公式可以較為準確地預測氣固流化床啟動階段內構件的受力特性。

圖4 流化床啟動階段內構件所受載荷實驗結果與CFD-DEM 模擬結果比較Fig.4 Comparison of the experimental and simulation results of the stress exerted on the baffle during start-up of the fluidized bed

4.2 顆粒碰撞模型參數對受力特性的影響

考慮到內構件表面載荷主要由固相顆粒與內構件表面的碰撞作用貢獻，在CFD-DEM 數值模擬方法中，離散顆粒與擋板壁面的碰撞作用主要由DEM 統計這些微觀物理量場信息。DEM 方法中，顆粒碰撞模型參數為碰撞作用的主要影響因素，顆粒楊氏模量、碰撞恢復系數、滑動摩擦因數和滾動摩擦因數這些碰撞參數對內構件表面與顆粒間的碰撞作用密切相關，并對擋板內構件表面載荷的影響尤為顯著。因此，考察固相顆粒的楊氏模量、摩擦因數和恢復系數對擋板表面受碰撞力的影響有重要意義。

顆粒楊氏模量是固相顆粒物性的重要參數，表征了顆粒表面的硬度。楊氏模量的大小取決于固相顆粒自身的性質，其數值越大顆粒越不容易發生變形，反映了固相顆粒材料的剛性[45]。此外，在計算機數值模擬方法中，離散顆粒的楊氏模量作為模型碰撞參數是CFD-DEM 數值模擬計算的一個重要參數。顆粒楊氏模量很大程度上影響數值模擬的計算量及計算效率，其決定了DEM 中計算的時間步長。在研究氣固兩相間的時均流動特性時，顆粒的楊氏模量對統計計算的流動參數無顯著影響，可以將楊氏模量數值設置小些以減少計算量，提高計算效率[46]。

對于密相床層的啟動階段，體系內顆粒-擋板壁面間的碰撞載荷統計則與楊氏模量的大小密切相關?；诒狙芯克x用的Hooke 模型，法向彈性系數kn與楊氏模量直接相關。因此，為了探究內構件的受力特性，需要研究楊氏模量的大小對氣固流化床啟動階段密相床層中擋板內構件表面受到的碰撞載荷的影響。固相顆粒物料仍為相同的B類石英砂顆粒，其楊氏模量的取值范圍介于0.1 GPa 和10 GPa 量級之間，其中前者為研究氣固流動特性統計時均參數常選用的經驗值，而后者是真實石英砂顆粒的楊氏模量量級，而擋板內構件的楊氏模量與顆粒的楊氏模量設置相同。本小節選用的顆粒楊氏模量具體數值及對應DEM 與CFD 中的時間步如表2 所示，模擬計算體系及其余模擬參數設置與前面小節模擬工況的設置一致。

表2 擋板床數值模擬參數Table 2 Simulation parameters for the baffled fluidized bed

圖5顯示了固相顆粒選取不同楊氏模量下密相床層中擋板內構件表面碰撞載荷隨時間的演化情況。發現流化床自啟動開始至出現第一個完整峰值脈沖載荷階段，因顆粒楊氏模量的不同，峰值載荷的大小與峰值出現的時間均受到影響發生相應變化，而后載荷的變化不再受顆粒楊氏模量的影響。對圖中啟動階段出現的第一個載荷峰值位置局部放大，如圖5 內附圖所示。（1）在顆粒初始化自由堆積成密相床層的階段(0 s 之前)，顆粒堆積過程中體系逐漸密實，新力鏈不斷生成。當床層高度穩定，此時體系足夠密實，力鏈充分發展幾乎不發生變化，擋板內構件表面載荷幾乎穩定不變[47]。顆粒堆積的隨機性與堆積時間的差異使得未通入氣體前擋板上、下表面受到的顆粒的擠壓碰撞作用不同，致使擋板內構件表面載荷存在一定的差異性。但是大體上，較大的顆粒楊氏模量使得顆粒-壁面間的法向接觸力增大，致使固定床狀態時擋板受到的碰撞作用載荷更大。（2）當體系內通入流化空氣，啟動階段開始至第一個峰值脈沖出現(0～0.12 s)，峰值出現的時間隨楊氏模量的增加而縮減，峰值載荷的大小隨顆粒楊氏模量的增大而減小，且當楊氏模量增大至10 GPa 量級時，此影響不再顯著，峰值大小不再發生明顯變化。說明顆粒楊氏模量的增加，加快了由固定床轉化為流化床的速度，楊氏模量大的顆粒床層內擋板表面載荷率先達到峰值。但由于楊氏模量增大，顆粒的表面硬度更大、變形程度減小，致使顆粒的碰撞重疊量減小，碰撞接觸時間明顯縮短，使得內構件表面的碰撞作用載荷減小。此外，楊氏模量大的顆粒床層因顆粒硬度大、變形小，在流體通入時會更容易打破“自鎖”狀態，固相顆粒較快發生松動而流化，可能也是使得擋板表面受碰撞載荷峰值較小的原因。（3）當擋板表面受到的第一個峰值脈沖結束后(0.12 s 之后)，內構件表面載荷變化基本不受顆粒楊氏模量的影響。體系進入正常流化階段后，顆粒楊氏模量大的床層中顆粒-壁面的碰撞作用更劇烈，發生碰撞更為頻繁，以及在彈性系數增加的共同影響下，在一定程度上抵消了因碰撞重疊量減小而使得內構件表面載荷減小的趨勢。

圖5 不同顆粒楊氏模量下顆粒床層內擋板受到載荷隨時間的變化Fig.5 Effect of particle Young's modulus on the stress exerted on the baffle immersed in the fluidized bed

流化床啟動階段，未流化密相床層會被氣節推起一定高度，使得擋板內構件受到較大的碰撞作用。而在這一短暫過程中，顆粒間的碰撞作用較為強烈，碰撞重疊量會更大，為保證模擬計算的精度，統計床層啟動階段顆粒最大重疊量與顆粒粒徑的比值隨時間的變化情況，如圖6 所示。從圖中可以看出，每個工況下顆粒最大重疊量與粒徑的比值都小于2%，說明啟動階段模擬計算的準確性。此外，隨著顆粒楊氏模量的增大，最大重疊量與粒徑的比值減小，也進一步說明了顆粒硬度增大，難以發生變形的特性。

圖6 不同楊氏模量下顆粒最大重疊量與粒徑比值隨時間的變化Fig.6 Variation of the ratio of the maximum particle-particle overlap to particle size with time under different Young's modulus

顆粒碰撞恢復系數是顆粒的重要物性參數之一，可以反映顆粒碰撞后恢復到變形前初始狀態的能力，即可以反映顆粒在接觸碰撞過程中能量耗散的強弱[48]。顆粒摩擦因數是固相顆粒物料的一個基本物性參數，與接觸面的粗糙程度有關。顆粒-顆粒與顆粒-壁面之間的摩擦會消耗顆粒的動能，進而影響體系內氣固相間的交換量以及固相顆粒在切向方向上的運動情況[49]。顆粒滾動摩擦因數作為固相顆粒物料的重要物理力學參數之一，也會影響顆粒體系內形成的微觀力學結構。顆粒間的滾動摩擦因數主要影響顆粒系統的堆積形態[50]，而在DEM 碰撞模型中，顆粒滾動摩擦因數的取值一般遠小于固相顆粒的直徑[51]。因此，依據單一變量原則，分別選取0.80、0.85、0.90 這三個值作為模擬所需的顆粒碰撞恢復系數；選取顆粒的滑動摩擦因數為0.30、0.35 和0.40；取顆粒的滾動摩擦因數為0.01、0.03 和0.05。其余參數的設置與之前保持一致，模擬結果如圖7所示，可以看出碰撞恢復系數、滑動摩擦因數和滾動摩擦因數對擋板內構件的受力影響并不顯著。

圖7 不同顆粒碰撞恢復系數(a)、滑動摩擦因數(b)、滾動摩擦因數(c)下顆粒床層內擋板受到載荷隨時間的變化Fig.7 Effect of particle restitution coefficient(a),particle sliding friction coefficient(b),and particle rolling friction coefficient(c)on the stress exerted on the baffle immersed in the fluidized bed

4.3 表觀氣速對受力特性的影響

流化床反應器的流化風速是控制流化床運行及氣固兩相流動的重要操作條件。因此，研究表觀氣速對擋板內構件受力特性的影響具有重要意義?；?.2 節關于碰撞模型參數對受力特性的影響研究，現選定一組關于顆粒屬性的最佳模型參數，即顆粒的楊氏模量為1×1010Pa，DEM 中時間步長為1×10-6s，CFD 中時間步長為1×10-5s，碰撞恢復系數為0.90，滑動摩擦因數為0.30，滾動摩擦因數為0.01。根據固相顆粒的屬性及最小流化氣速，現選取用于模擬研究的表觀氣速的一組值0.5、0.6、0.9 m/s。

圖8 為不同表觀氣速下，流化床啟動階段密相床層中擋板內構件受到的碰撞載荷瞬態變化情況。模擬結果顯示，對于初始堆積床層，表觀氣速增大，相當于給予流化床的氣速增量增加，擋板內構件的受力載荷強度也相應增大。當表觀氣速增加的幅度越大，內構件表面載荷的數值增長的幅度也越大。出現這一現象的原因在于表觀氣速增大，氣速增量越大，流化床底部形成的氣節高度更高，床層中固相顆粒的動量越大，對擋板內構件表面的撞擊更為劇烈，因此載荷強度在數值上表現出更大的情況。說明在流化床床層啟動階段，表觀氣速的大小是影響擋板內構件受力的主要因素之一。

圖8 不同表觀氣速下顆粒床層內擋板受到載荷隨時間的變化Fig.8 Effect of superficial gas velocity on the stress exerted on the baffle immersed in the fluidized bed

上述模擬結果說明在實際流化床啟動操作過程中，設置較小的表觀氣速可明顯減小內構件受到的碰撞力。因此，要避免在較大的氣速通入下完成床層啟動操作，防止擋板內構件因瞬時受到較大的作用力而發生強度破壞，導致反應器裝置無法正常運行。

4.4 顆粒粒徑對受力特性的影響

流化床體系內，固相顆粒粒徑的改變可使系統的最小流化速度發生變化，從而影響體系內的流態轉變及顆粒碰撞動力學[49]。因此，采用不同粒徑的顆粒物料作為床料進行數值模擬，基于單一變量原則，在表觀氣速分別按統一設定的工況和以最小流化氣速倍數設定的工況下，研究顆粒粒徑變化對密相床層中擋板內構件受力特性的影響。模擬研究所選用的固相顆粒僅粒徑不同，其余的顆粒物性參數均一致（如顆粒數目）。

基于單一變量原則，將流化床內的流化氣速采用統一設定，并將其設置為0.6 m/s，以此來研究相同表觀氣速下顆粒尺寸對內構件受力特性的影響。本節關于顆粒屬性的模型參數設置與4.3 節一致，其余參數參考表1設置。

圖9 顯示了表觀氣速統一設定為0.6 m/s 工況下，不同粒徑的顆粒床層啟動階段密相床層中擋板內構件表面載荷強度的瞬態演變情況。從圖中可看出隨著顆粒粒徑的增大，內構件表面受的碰撞載荷強度整體有減小趨勢，尤其在第一個完整峰值載荷出現的區間，載荷強度有明顯減小的趨勢，在0.06 s 位置表現更為明顯。出現這一現象的原因可能在于，隨著顆粒尺寸的增加，體系內的最小流化氣速相應增加，而在相同表觀氣速下粒徑較大的固相顆粒床層易出現流化不良的現象，床層中氣泡的尺寸相對較小，氣固兩相的運動及相間作用程度減弱，擋板內構件受到的載荷作用也相應減小，這與前期實驗中的現象基本一致。

5 結 論

本研究基于CFD-DEM 數值模擬方法，基于劉對平[4]前期實驗所用的三維方形冷模流化床裝置作為模擬體系，以工業上廣泛應用的單個斜片擋板內構件作為研究對象，建立了擋板內構件表面受力載荷強度的統計計算方法，并借助此統計計算公式分析研究了氣固物性對擋板床啟動階段密相床層中單個水平斜片擋板內構件的受力特性的影響規律。主要結論如下。

（1）基于CFD-DEM 模擬方法建立了氣固流化床中浸沒在密相床層的擋板內構件表面受顆粒碰撞載荷分布的統計計算公式，可用于定量分析內構件表面受到的動態碰撞載荷信號。通過與課題組前期實驗測量數據對比，建立的內構件表面載荷統計計算公式可以定性和(半)定量復現實驗中流化床啟動階段內構件受力特性的變化，為研究內構件在床層啟動階段流型轉變過程中出現這一非正常受力現象的內在機理奠定基礎。

（2）顆粒碰撞模型參數中，楊氏模量是影響內構件受力特性的關鍵性因素，而碰撞恢復系數、滑動摩擦因數和滾動摩擦因數對擋板內構件的受力影響并不顯著。研究發現流化床自啟動開始至出現第一個完整峰值脈沖載荷階段，不同的顆粒楊氏模量下，內構件受力的峰值載荷大小與峰值出現的時間均受到影響發生相應變化（隨著楊氏模量的增大，峰值載荷先逐漸變小后趨于穩定），而后階段載荷的變化幾乎不再受顆粒楊氏模量的影響。此外，在受楊氏模量影響階段，隨著楊氏模量的增大，第一個完整峰值載荷相應減小，當楊氏模量增大至顆粒物理真實值，峰值載荷不再發生變化。

（3）氣固流化床體系內表觀氣速是影響內構件受力的主要因素之一。表觀氣速越大，密相床層中的擋板內構件受力載荷強度越大，與實驗結論一致。為保障工業上反應器裝置的正常運行和內部構件及其支撐結構的穩定性，要在相對較小的氣速增幅下完成流化床的啟動操作過程。