壓縮空氣泡沫管路輸送的模擬計算

林全生，高 紅，宋 波

（1. 天津大學 化工學院，天津 300350；2. 公安部天津消防研究所，天津 300381）

［關鍵字］泡沫；數值模擬；黏度；沿程阻力系數；Fluent軟件

壓縮空氣泡沫滅火劑是一類主要用于撲救固體可燃物火災的新型高效泡沫滅火劑，具備產生高穩定泡沫的性能，可用于石油化工可燃性液體危險品儲罐火災時的比鄰儲罐防護［1-4］。該類滅火劑在固定式壓縮空氣泡沫系統中具有突出的滅火效果。但泡沫在管道運輸過程中，以水為基礎的現有理論難以滿足要求，需要找到適合壓縮空氣泡沫流體流動的計算方法，為工程設計提供準確依據。針對泡沫，利用單相流對其進行研究的較多。葛曉霞等［5］通過搭建實驗平臺來觀察泡沫在管內的流態，給出了壓縮空氣泡沫流型的類型，并通過理論分析得出了模型參數的算法，確定了摩阻系數及壓力損失的計算方程。Sun等［6］通過實驗發現泡沫流體是剪切稀化的冪率流體，并得到了摩阻系數與廣義雷諾數的關系。陳旸等［7］開展了壓縮空氣泡沫管網輸送實驗，研究了泡沫流量、氣液比、泡沫種類對輸送管路壓力變化的影響及摩擦阻力損失。Hu等［8］對氣液兩相流在90°彎管中的流動情況進行了模擬，得到了一定氣液比下雙相流彎管壓力分布規律。Supa-Amornkul等［9］分別采用單相流和雙相流對彎管壓降進行數值模擬，氣體體積分數在50%以上的情況下，沿長度方向單相流壓力分布與實驗結果差別不大，相比之下，多相流并不能很好地預測彎管壓降。林全生等［10］利用上海消防總隊的實驗數據，對泡沫在不同管徑、不同長度的管路中的流動阻力損失進行了數值模擬，并得出了范寧摩擦因子與廣義雷諾數的擬合關系式。Fluent數值模擬已在工程上有很多應用［11-14］，它在消防泡沫輸運方面的研究還有很大的探索空間。

本工作利用Fluent數值模擬計算軟件對壓縮空氣泡沫在管道內的流動進行了數值模擬研究，考慮了發泡倍數和流量的影響，并與相關實驗數據進行了對比，通過擬合得出沿程阻力系數的顯性表達式，為壓縮空氣泡沫的工程計算提供指導。

1 數值模擬計算方法

1.1 幾何模型和網格劃分

基于有限元方法，根據文獻［7］，管路設備尺寸如圖1所示，利用Solidworks軟件建立三維幾何模型。利用ANSYS Workbench16.0中的mesh模塊對幾何模型進行四面體網格劃分，網格數量為113萬，并對幾何模型的壁面設置膨脹層進行邊界層加密處理，如圖2所示。

圖1 實體管路示意圖［7］Fig.1 Diagram of solid piping［7］.

圖2 幾何模型網格劃分Fig.2 Mesh division of geometric model.

1.2 計算模型假設

葛曉霞等［5］對壓縮空氣泡沫進行了流型的劃分，將泡沫在管內的流動形態按含氣率從低到高的順序依次劃分為泡狀流、塞狀流、層狀流、波狀流、彈狀流及環狀流。Rouhani等［15］對兩相流的流型進行了劃分，認為水平流動主要包括平滑流、波浪流、柱塞流和彌散流等。消防領域應用的泡沫氣體體積分數一般大于80%，泡沫流動過程中是液膜包圍著氣體運動，泡沫形態與Rouhani等提出來的彌散流相似，泡沫在管道內的流動是充分發展的，且管道內的泡沫湍流程度較高，導致析出的液體重新混合，從而使泡沫一直保持一個較為均勻的狀態，不會發生兩相分離現象［10］。故而，假設泡沫在同一管道截面為均勻流動，進一步將泡沫簡化為均質單相流體［10，16］。實際工況下，對于泡沫在管路中的流動不研究熱量的傳遞，因此模擬過程中不計流動過程中的能量損失。假設不存在壁面滑移現象。

1.3 模型的選擇及邊界條件的設定

湍流在實際工業中的存在較為普遍，實際工程應用中，泡沫流速較大，大部分的泡沫流體都處于湍流狀態。Spalart-Allmaras（S-A）模型適用于雷諾數較低的湍流情況，因此采用S-A單方程模型模擬壓縮空氣泡沫在管道內的流動。此外網格需滿足黏性影響區流體流動求解要求，故壁面處網格需要局部加密［10］。選用SIMPLE算法，針對文獻［7］中1%的A類泡沫氣液比7∶1和10∶1的兩種發泡倍數，在不同流量下的圓形管路中的流動進行模擬研究。本工作采用體積流量換算出的速度作為入口邊界條件，氣液比為7∶1時入口流量分別為846，1 660，2 420，3 300，3 930，5 060 L/min，氣液比為10∶1時入口流量分別為1 240，2 230，3 300，4 420 L/min。出口邊界條件為壓力出口，由于文獻［7］中管路出口有球形閥，球形閥出口邊界條件的壓力值為標準大氣壓。壁面邊界條件為壁面無滑移條件。

1.4 控制方程

本工作不考慮熱量的損失，不考慮能量方程。采用S-A單方程模型求解泡沫在管路中的流動問題?；究刂品匠贪ㄟB續性方程、動量方程以及S-A方程。

連續性方程見式（1）:

式中，ρ為流體密度，kg/m3；t為時間，s；Vx為x方向的流體流動速度，m/s；Vy為y方向的流體流動速度，m/s。

動量方程見式（2）:

式中，p為流體流動中靜壓，Pa；ui和uj為速度張量，m/s；τij為應力張量，Pa；ρgi與Fi分別是i方向上重力體積力和外部體積力，N；xi和xj為坐標張量。

上式中應力張量形式可表達為式（3）:

式中，τij是由于流體中分子黏性而在微元體表面產生的黏性應力τ的分量；μa為流體黏性系數；δij為克羅內克δ符號。

圖3 氣液比為7∶1時文獻［7］實驗值與模擬值結果的對比Fig.3 The comparison of experimental value in literature［7］ and simulated result when the gas-liquid ratio is 7∶1.

1.5 Fluent數值模擬

根據文獻［7］中的實驗數據，需要得出模擬過程中第1個壓力傳感器與第6個壓力傳感器處的壓力值。利用Fluent軟件在傳感器1和傳感器6處建立截面，提取出該位置的壓力值。利用數據點1與數據點6的壓力值的差值作為阻力損失，并分別將數據點1和數據點6的數值模擬值與實驗值進行比較。

2 結果與討論

2.1 壓縮空氣泡沫壓力損失模擬分析

將文獻［7］中阻力損失數據與模擬結果進行對比。氣液比分別為7∶1和10∶1時的壓力對比及阻力損失誤差見圖3和圖4，泡沫流量分別為3 000 L/min和5 000 L/min時的壓力對比及阻力損失誤差見圖5和圖6。由圖3～6可看出，模擬結果的相對誤差基本在10%以內，數值模擬的結果在工程誤差允許的范圍內。

圖4 氣液比為10∶1時文獻［7］實驗值與模擬值結果的對比Fig.4 The comparison of experimental value in literature［7］ and simulated result when the gas-liquid ratio is 10∶1.

圖5 流量為3 000 L/min時不同發泡倍數下文獻［7］實驗值與模擬值結果的對比Fig.5 The comparison of experimental value in literature［7］ and simulated result at different foaming multiples with flow rate of 3 000 L/min.

圖6 流量為5 000 L/min時不同發泡倍數下文獻［7］實驗值與模擬值結果對比Fig.6 The comparison of experimental value in literature［7］ and simulated result at different foaming multiples with flow rate of 5 000 L/min.

2.2 數值模擬過程的黏度分析

由文獻［7］中表觀黏度隨剪切速率的變化關系可以得知，由于實驗條件的限制只能測出剪切速率230 s-1以下的泡沫黏度，然而實際工程中泡沫流體速度較大，泡沫流體的剪切速率也較大。通過計算可以得出泡沫流體在較高流速下的剪切速率，具體方法如下。

非牛頓流體在管內流動的基本假設:流體已充分發展，管壁處不發生滑移；不考慮流變性與時間的關系；流體在管內等溫流動，且流動邊界直徑不變。根據以上假設，有非牛頓流體管內流動公式，即Rabinowitsch-Mooney方程，見式（5）:

式中，γw為管壁處的剪切速率，s-1；τw為管壁處的剪切應力，Pa；u為流體的流動速度，m/s；D為管路內徑，m；n為流變指數。

由Rabinowitsch-Mooney方程可知，只要確定了流變指數n的值，就可以把流體速度和剪切速率聯系起來。將文獻［7］中的流量轉換成速度u，通過對文獻［7］中黏度和剪切速率的關系進行擬合得出流變指數n，由式（5）即可以計算出實驗條件下的剪切速率。根據陳旸等［7］的研究，壓縮空氣A類泡沫的流變規律符合Herschel-Bulkley模型，見式（6）:

公式兩邊同時除以剪切速率γ，得式（7）:

式中，K為流變系數；τ為剪切應力，Pa；τ0為剪切速率為零時的剪切應力，Pa；γ為剪切速率，s-1。從而計算得出高剪切速率下的泡沫流體表觀黏度。

泡沫流體是非牛頓流體，具有剪切變稀的性質。通過對不同氣液比及不同流量下的泡沫在直徑不變的管路中流動進行數值模擬，最終誤差在10%之后再確定該條件下的黏度，從而得出模擬條件下不同氣液比及不同流量下壓縮空氣泡沫的黏度。將計算所得泡沫流體的表觀黏度與模擬條件得出的表觀黏度進行對比，結果見圖7。

圖7 計算值與模擬值的對比Fig.7 Comparison of calculated and simulated values.

由圖7可看出，模擬與計算所得表觀黏度值吻合相對較好。因此Rabinowitsch-Mooney方程能夠較準確地預測高流速下的泡沫剪切速率。

2.3 沿程阻力系數分析

根據1.3節中的邊界條件，對泡沫在管路中的流動進行數值模擬。水平管路總的壓力損失Δp可以表示為式（8）:

式中，Δpf為沿程阻力損失，Pa；Δpa為流體流動時加速度壓力降，Pa。

當流體流動是充分發展的穩態流動時，流體的加速度壓力降是可以忽略不計的。將式（8）與達西公式結合可以得出沿程阻力系數與壓力損失的關系，見式（9）:

式中，λ為沿程阻力系數；L為管路長度，m。

從式（9）可知，沿程阻力系數、管道長徑比以及泡沫密度、流速都是影響流動壓降的重要因素。Sun等［6］總結了沿程阻力系數的關系式，見式（10）:

式中，T為溫度，℃。

在不考慮溫度的情況下，沿程阻力系數可表示為雷諾數的函數，即λ=f（Re），如圖8所示。

對計算結果進行擬合，得到式（11）。

對于直徑為80 mm的管路，氣液比為7∶1和10∶1的不同發泡倍數的泡沫沿程阻力系數與雷諾數的關系，它們的理論計算值與數值模擬值具有較好的一致性。林全生等［10］對泡沫流量不變、發泡倍數不變，改變管路直徑的情況進行了研究，得出了沿程阻力系數與廣義雷諾數的關系，但是沒有考慮發泡倍數與流量的影響，因此僅適用于發泡倍數為6.7倍的情況，對于泡沫發泡倍數變化的情況則不適用。根據文獻［7］中的數據，對不同發泡倍數的泡沫在直徑不變的管路中的流動進行了研究，所以得出的沿程阻力系數與廣義雷諾數的關系僅適用于管徑不變、泡沫發泡倍數與泡沫流量變化的情況。

圖8 沿程阻力系數與雷諾數的關系Fig.8 The relation between drag coefficient along the path(λ) and Re.

2.4 不同流量及不同發泡倍數下出口壓力的變化

氣液比為7∶1和10∶1時泡沫出口壓力隨流量的變化如圖9所示。

圖9 不同氣液比下出口壓力隨流量的變化Fig.9 The variation relationship between outlet pressure and flow rate under different gas-liquid ratios.

從圖9可發現，在實驗所選用的流量下，出口壓力隨流量的增加而增大，這是因為泡沫流量增加，泡沫流體對出口的沖擊力增大，泡沫出口壓力增大。泡沫出口壓力與測試點處的壓力具有相同的變化趨勢。這對實際工程應用具有一定的指導意義。

3 結論

1）數值模擬所得管路中的壓力損失與實驗值具有較好的一致性，模擬結果誤差在10%以內。

2）根據Rabinowitsch-Mooney方程，可以計算出高流速下的剪切速率；再根據黏度隨剪切速率的變化關系，計算得出高流速下的泡沫黏度。根據Rabinowitsch-Mooney方程計算出的泡沫黏度與數值模擬所得泡沫黏度吻合較好。

3）利用達西公式分析并建立了沿程阻力系數與雷諾數的關系式λ= 2.998 3/Re0.5909，氣液比7∶1和10∶1的不同發泡倍數的泡沫沿程阻力系數與雷諾數的關系的文獻實驗值與數值模擬值具有較好的一致性。

4）模擬結果表明，管路出口處的壓力隨流量的增加而增加，并且在不同測壓點具有相同的變化趨勢，這對實際工程應用具有一定的指導意義。

符 號 說 明

CS-A模型中的常數，C=0.622

D管路內徑，m

Fi外部體積力，N

G湍流黏性產生項，kg/（m3·s）

K流變系數

L管路長度，m

n流變指數

p流體流動中靜壓，Pa

Δp水平管路總的壓力損失，Pa

Δpa加速度壓力降，Pa

Δpf沿程壓力降，Pa

Re雷諾數

T溫度，℃

t時間，s

u流體的流動速度，m/s

ui速度張量，m/s

uj速度張量，m/s

vxx方向的流體流動速度，m/s

vyy方向的流體流動速度，m/s

xi坐標張量

xj坐標張量

Y壁面阻擋與黏性阻尼引起的湍流黏性的減少項，kg/（m3·s）

γ剪切速率，s-1

γw管壁處的剪切速率，s-1

δij克羅內克δ符號

λ沿程阻力系數

μ流體黏度，Pa·s

μa流體黏性系數

ρ流體密度，kg/m3

ρgi重力體積力，N

σS-A湍流模型中的常數，σ=2/3

τ剪切應力，Pa

τ0剪切速率為零時的剪切應力，Pa

τw管壁處的剪切應力，Pa

τij應力張量，Pa

參 考 文 獻

［1］ 陸強，包志明，陳濤，等. A類泡沫滅火劑滅火性能試驗研究［J］. 消防科學與技術，2013，32（2）:177-179.

［2］ 初迎霞，王海明，喬啟宇，等. 壓縮空氣泡沫系統在林火撲救中的應用［J］. 林業機械與木工設備，2004，32（8）:33-36.

［3］ 傅學成，葉宏烈，包志明，等. A類泡沫滅火劑的發展與瞻望［J］. 消防科學與技術，2008，27（8）:590-592.

［4］ 傅學成，包志明，陳濤，等. 壓縮空氣泡沫的隔熱防護性能研究［J］. 消防科學與技術，2009，28（3）:204-207.

［5］ 葛曉霞，高健，趙昊. 壓縮空氣泡沫管內流動特性分析［J］.消防科學與技術，2016，35（10）:1408-1410.

［6］ Sun Xiao，Liang Xiaobing，Wang Shuzhong，et al. Experimental study on the rheology of CO2viscoelastic surfactant foam fracturing fluid［J］. J Pet Sci Eng，2014，119（3）:104-111.

［7］ 陳旸，陳濤，胡成，等. 壓縮空氣泡沫傳輸阻力損失研究［J］. 消防科學與技術，2017，36（10）:1418-1420.

［8］ Hu Dongfang，Huang Zhengliang，Sun Jingyuan，et al. Numerical simulation of gas-liquid flow through a 90° duct bend with a gradual contraction pipe［J］. J Zhejiang Univ Sci A，2017，18（3）:212-224.

［9］ Supa-Amornkul S，Steward F R，Lister D H. Modeling twophase flow in pipe bends［J］. J Pressure Vessel Technol，2005，127（2）:204-209.

［10］ 林全生，張猛，宋波，等. 壓縮空氣A類泡沫在水平管道內流動研究［J］. 消防科學與技術，2017，36（9）:1265-1268.

［11］ 張仂，谷芳. 基于計算流體力學熱流固耦合仿真的換熱器折流板結構優化［J］. 石油化工，2012，41（6）:682-687.

［12］ Pravakar M，Khakhar D V. Hydraulic resistance of rigid polyurethane foams. I. Effect of different surfactants on foam structure and properties［J］. J Appl Polym Sci，2004，93（6）:2821-2829.

［13］ 甘德清，高鋒，陳超，等. 管道輸送高濃度全尾砂充填料漿的阻力損失研究［J］. 礦業研究與開發，2016，36（1）:94-98.

［14］ 王新民，丁德強，吳亞斌，等. 膏體充填管道輸送數值模擬與分析［J］. 中國礦業，2006，15（7）:57-59.

［15］ Rouhani S Z，Sohal M S. Two-phase flow patterns:A review of research results［J］. Prog Nucl Energy，1983，11（3）:219-259.

［16］ 王勇凱，高紅，宋波，等. 壓縮空氣 A 類泡沫水平管路壓降試驗及數值模擬［J］. 化工學報，2018，69（10）:4184-4193.