孔隙尺度下超臨界CO2驅(qū)水兩相流數(shù)值模擬

姜水生，趙萬東，張瑩，李培生，王昭太，鐘源

（南昌大學(xué)機(jī)電工程學(xué)院，江西 南昌 330031）

孔隙尺度下超臨界CO2驅(qū)水兩相流數(shù)值模擬

姜水生，趙萬東，張瑩，李培生，王昭太，鐘源

（南昌大學(xué)機(jī)電工程學(xué)院，江西 南昌 330031）

以咸水層封存二氧化碳（CO2）為研究背景，利用流體體積函數(shù)（Volume of Fluid，VOF）多相流模型建立孔隙尺度多孔介質(zhì)計(jì)算模型，研究了超臨界二氧化碳 （Sc-CO2）注入到含水多孔介質(zhì)中的兩相運(yùn)移規(guī)律。對(duì)比分析了毛細(xì)管數(shù)、地質(zhì)封存壓力、Sc-CO2注射溫度、兩相表面張力系數(shù)、接觸角等因素對(duì)兩相運(yùn)移速率以及驅(qū)替效率的影響，同時(shí)將不同毛細(xì)管數(shù)下的驅(qū)替效率與實(shí)驗(yàn)數(shù)據(jù)進(jìn)行了對(duì)比。研究結(jié)果表明，隨著毛細(xì)管數(shù)的增大，驅(qū)替效率先減小然后趨于穩(wěn)定，數(shù)值模擬與實(shí)驗(yàn)數(shù)據(jù)吻合良好；在同一孔隙率下，在壁面表現(xiàn)為親水性時(shí)，壁面潤(rùn)濕性越好，驅(qū)替速率越快，但驅(qū)替效率有所下降；同時(shí)毛細(xì)管數(shù)越小、地質(zhì)封存壓力越低、注射溫度越高、張力系數(shù)越小驅(qū)替速率越快，且驅(qū)替效率越高。

超臨界二氧化碳；咸水層封存；孔隙尺度；數(shù)值模擬；兩相流；毛細(xì)管數(shù)；接觸角

隨著溫室效應(yīng)問題的日益嚴(yán)重，溫室氣體的排放與控制成為當(dāng)今學(xué)者們研究的熱點(diǎn)。碳元素的捕集與封存（carbon capture and storage，簡(jiǎn)稱CCS）[1-2]是減少溫室氣體的重要措施之一，而CCS中地質(zhì)封存被認(rèn)為是中短期內(nèi)減少大氣中CO2含量的最有效方法之一。地質(zhì)封存分為深部咸水層封存（簡(jiǎn)稱鹽水層封存）[3-4]、枯竭油氣藏封存[5]以及不可開采的煤層封存，CO2封存示意圖如圖1所示。其中咸水層封存因?yàn)楹叵滤约翱蓾B透的多孔介質(zhì)，且廣泛分布于世界各地，故被認(rèn)為是最具潛力的封存地點(diǎn)。超臨界二氧化碳（Sc-CO2）擴(kuò)散性好且動(dòng)力黏度小，故將Sc-CO2注入到咸水層中成為封存CO2研究的熱點(diǎn)。

圖1 CO2地質(zhì)封存示意圖

近年來國(guó)內(nèi)外學(xué)者對(duì)Sc-CO2咸水層封存運(yùn)移機(jī)理進(jìn)行了不同程度的研究，SUEKANE等[6]將Sc-CO2注入到使用玻璃珠堆積形成的人工多孔介質(zhì)中，使用核磁共振（MRI）技術(shù)進(jìn)行了可視化實(shí)驗(yàn)研究；WANG等[7]通過實(shí)驗(yàn)，使用連續(xù)和不連續(xù)的注射速率將Sc-CO2注射到二維均質(zhì)咸水多孔介質(zhì)中，研究了毛細(xì)管數(shù)對(duì)驅(qū)替過程中Sc-CO2飽和度的影響；TREVISAN等[8-9]通過實(shí)驗(yàn)構(gòu)建均質(zhì)和非均質(zhì)多孔介質(zhì)，進(jìn)行巖心驅(qū)替實(shí)驗(yàn)，并使用X射線衰減分析流域飽和度，確定毛細(xì)阻力的影響因素；SASAKI等[10]通過直接數(shù)值模擬方法對(duì)不同壓力下的Sc-CO2注入到地下巖石進(jìn)行了研究；CHOI等[11]為了研究CO2封存過程中的安全性，建立了孔隙尺度下的多孔介質(zhì)模型，使用CFD軟件STAR-CCM研究了不同注射溫度下的CO2驅(qū)替過程中的流動(dòng)和傳熱過程。姜培學(xué)課題組[12]采用核磁共振設(shè)備，通過實(shí)驗(yàn)研究了Sc-CO2和鹽水流過玻璃珠過程中浮升力和孔隙度的影響，同時(shí)使用可視化技術(shù)研究了巖心中超臨界壓力下CO2-水相對(duì)滲透率與飽和度的關(guān)系；XU等[13]在孔隙尺度條件下直接求解N-S方程方法對(duì)比研究了Sc-CO2驅(qū)水和空氣驅(qū)水過程，但只模擬了一種方案下的驅(qū)替過程；羅庶等[14]使用數(shù)值模擬方法研究了地質(zhì)封存Sc-CO2突破壓力梯度的問題，但對(duì)多種因素影響下的驅(qū)替過程中兩相界面變化、突破過程、Sc-CO2飽和度等研究較少；蔣蘭蘭[15]通過MRI技術(shù)以及CT成像技術(shù)，對(duì)玻璃砂多孔介質(zhì)進(jìn)行了單相和多相驅(qū)替實(shí)驗(yàn)，并使用瞬態(tài)和穩(wěn)態(tài)法研究了相對(duì)滲透率飽和曲線。

國(guó)內(nèi)外學(xué)者對(duì)Sc-CO2驅(qū)水過程進(jìn)行了大量的實(shí)驗(yàn)研究，然而在實(shí)驗(yàn)過程中，MRI在壓力的選定上存在限制，且實(shí)驗(yàn)周期長(zhǎng)，成本大；在驅(qū)替實(shí)驗(yàn)中，由于其構(gòu)造的多孔介質(zhì)模型單一，表面張力以及接觸角等因素較難改變；同時(shí)由于出口末端效應(yīng)和入口段效應(yīng)的存在[15]，不能精確地測(cè)量多孔介質(zhì)兩端的壓力，因此存在一定的實(shí)驗(yàn)誤差。而數(shù)值模擬方法在選取不同的孔隙率、毛細(xì)管數(shù)、注射壓力、表面張力、接觸角等參數(shù)較容易。故本文使用流體體積函數(shù)（VOF）方法對(duì)Sc-CO2驅(qū)水進(jìn)行了孔隙尺度單相驅(qū)替研究，研究不同注射壓力和溫度等因素對(duì)兩相界面變化、突破過程、驅(qū)替速率、驅(qū)替效率以及殘余水飽和度的影響，將有助于CO2地質(zhì)封存場(chǎng)地的選取以及Sc-CO2注射條件的選取。

1 物理模型、控制方程和數(shù)值求解方法

1.1 物理模型

描述真實(shí)多孔介質(zhì)三維模型十分困難，而通過孔隙尺度下的有限微小單元可以描述不同孔隙率的多孔介質(zhì)；假設(shè)固體顆粒大小均勻，同時(shí)顆粒排布均勻，如實(shí)驗(yàn)中人工填充的玻璃珠模型。改變兩個(gè)簡(jiǎn)單的不同直徑大小的圓球來實(shí)現(xiàn)流體域的建立，例如簡(jiǎn)單立方體（simple cube）、體積中心立方體（body-centered cube）以及面中心立方體（face-centered cube）方式。XU和JIANG[16]已經(jīng)通過數(shù)值模擬研究了不同排列方式下，考慮壁面的剪切速度以及從實(shí)驗(yàn)方面對(duì)比3種模型的誤差，其中體積中心立方體模型的誤差小于5%。因此本文使用體積中心立方體模型來描述孔隙尺度模型，將孔隙結(jié)構(gòu)中固體域分離即可得出流體計(jì)算域，計(jì)算域孔隙尺度單元以及尺寸見圖2與表1。求解域離散成CFD計(jì)算網(wǎng)格，并將流體域分為質(zhì)量流量入口、壓力出口、壁面、對(duì)稱邊界。網(wǎng)格邊界層局部放大圖以及邊界條件見圖3，流體計(jì)算域離散網(wǎng)格使用ANSYS生成。

圖2 孔隙尺度結(jié)構(gòu)單元以及體積中心立方體模型

表1 單元尺寸以及孔隙度

圖3 離散網(wǎng)格邊界條件及橫截面局部網(wǎng)格放大圖

1.2 控制方程

Sc-CO2驅(qū)替含水多孔介質(zhì)是一個(gè)復(fù)雜的多相流過程，涉及到相與相之間的界面作用力，如表面張力的作用，同時(shí)還有Sc-CO2、水以及壁面之間的三相接觸角的問題。VOF方法已經(jīng)廣泛應(yīng)用于多相流仿真，如模擬潰壩、沸騰與冷凝、氣泡的運(yùn)動(dòng)、自由界面等模擬[17-18]。VOF基于追蹤技術(shù)應(yīng)用于固定的歐拉網(wǎng)格上的自由界面，求解模型中引入一個(gè)體積比函數(shù)α，假設(shè)αq為其中主相的體積分?jǐn)?shù)，也即該相占整個(gè)網(wǎng)格體積的體積分?jǐn)?shù)[19]。通過此方法能較好地實(shí)現(xiàn)對(duì)自由面以及兩相流內(nèi)部運(yùn)動(dòng)交界面的追蹤，因此可以實(shí)現(xiàn)Sc-CO2對(duì)驅(qū)替水過程中的內(nèi)部運(yùn)移界面的追蹤。

多相界面的追蹤通過求解連續(xù)性方程主相來完成，對(duì)于其中的某一相連續(xù)性方程表達(dá)如式(1)[19]。

式中，mpq為從p到q的質(zhì)量轉(zhuǎn)移，相反mpq是p到q的質(zhì)量轉(zhuǎn)移，同時(shí)首相的體積分?jǐn)?shù)滿足如下式(2)定值。

其中體積分?jǐn)?shù)方程通過顯示時(shí)間離散求解，在顯示求解過程中，通過有限差分插值格式應(yīng)用于體積分?jǐn)?shù)計(jì)算的前一個(gè)時(shí)間步，計(jì)算格式如式(3)。

式(3)中，n+1為當(dāng)前的計(jì)算時(shí)間步；n代表上一個(gè)時(shí)間計(jì)算步。

通過一個(gè)動(dòng)量方程求解整個(gè)計(jì)算域，得到的速度、溫度場(chǎng)應(yīng)用于整個(gè)流場(chǎng)，式(4)是動(dòng)量方程的表達(dá)式。

在控制體的輸運(yùn)方程中密度和黏度由混合相確定，在兩相求解模型中，密度和黏度如式(5)、式(6)定義。

求解過程中添加能量項(xiàng)、能量方程的計(jì)算如式(7)、式(8)。

式(7)中，ρ為式(5)計(jì)算得到的平均密度；keff為有效熱導(dǎo)率。

求解過程引入表面張力，也即是CSF模型，該模型把表面張力看作求解過程中的體積力添加在動(dòng)量方程中，引入的體積力如式(9)。

式(9)中，σqp為兩相表面張力系數(shù)；ρ為式(5)計(jì)算得到的平均密度。

三相界面接觸角采用固定接觸角，接觸角定義如式(10)。

式(10)中，θw為壁面接觸角；為分別為壁面法向和切向單位向量。

1.3 數(shù)值求解方法與邊界條件

由于多孔介質(zhì)結(jié)構(gòu)的復(fù)雜性，為了表征邊界層網(wǎng)格，因此離散網(wǎng)格采用四面體網(wǎng)格。使用Fluent16.0中的VOF方法模擬Sc-CO2驅(qū)替含水多孔介質(zhì)過程，求解模型設(shè)定為瞬態(tài)VOF模型，考慮體積力打開implicit body force選項(xiàng)；由于驅(qū)替過程雷諾數(shù)較小，因此求解模型選擇層流模型；模擬過程中有能量的傳遞，故加入能量方程；考慮到瞬態(tài)求解的收斂性，科拉數(shù)設(shè)定為0.5，壓力與速度的耦合采用PISO算法，動(dòng)量方程采用二階迎風(fēng)格式，同時(shí)采用幾何重構(gòu)的方法獲取Sc-CO2驅(qū)替的界面形狀；殘差收斂設(shè)定為1×10–4，求解過程中設(shè)定時(shí)間步長(zhǎng)為1×10–6s。入口邊界條件設(shè)定為Sc-CO2質(zhì)量流量入口，質(zhì)量流量為5×10–4g/s，出口設(shè)定為壓力出口，靜壓為0，壁面選擇無滑移邊界條件，環(huán)境變量根據(jù)選取的地質(zhì)封存深度來決定，含水多孔介質(zhì)初始溫度設(shè)定為333K，其他邊界設(shè)定為對(duì)稱邊界條件。

VOF材料的添加選取水和Sc-CO2。對(duì)于Sc-CO2材料的添加，可通過界面操作命令來實(shí)現(xiàn)，界面操作命名中輸入NIST Real Gas Models，選取單組分Sc-CO2材料，并在界面命名中創(chuàng)建計(jì)算節(jié)點(diǎn)表。本次數(shù)值模擬方案以及Sc-CO2、水的物性參數(shù)分別見表2和圖4、圖5。

1.4 數(shù)據(jù)處理方法

對(duì)于評(píng)定Sc-CO2驅(qū)替過程參數(shù)有滲透率、絕對(duì)滲透率、相對(duì)滲透率等。對(duì)于服從達(dá)西定律的多孔介質(zhì)流動(dòng)，滲透率可以計(jì)算出，而絕對(duì)滲透率是在多孔介質(zhì)內(nèi)的單相流動(dòng)評(píng)定參數(shù)。對(duì)多孔介質(zhì)內(nèi)的絕對(duì)滲透率定義如式(11)[15]。

表2 數(shù)值模擬方案及相應(yīng)參數(shù)的選取

圖4 Sc-CO2物性參數(shù)隨溫度和壓力變化關(guān)系

式中，Q為流量；L為長(zhǎng)度；μ為流體動(dòng)力黏度；ΔP為壓降。多孔介質(zhì)滲透率只與多孔介質(zhì)本身有關(guān)，而與流體屬性無關(guān)，本文采取穩(wěn)定流動(dòng)獲取多孔介質(zhì)滲透率，計(jì)算出的滲透率與實(shí)驗(yàn)數(shù)據(jù)吻合[15]。

計(jì)算過程中采用量綱為1毛細(xì)管數(shù)來量化不同方案，毛細(xì)管數(shù)定義如式(12)[15]。

式中，μ為Sc-CO2動(dòng)力黏度，Pa·s；V為速度，m/s；σ為兩相表面張力系數(shù)，N/m。

計(jì)算過程中水和Sc-CO2的飽和度定義如式(13)、式(14)。

2 結(jié)果與討論

2.1 驗(yàn)證性計(jì)算

2.1.1 網(wǎng)格無關(guān)性驗(yàn)證

數(shù)值模擬過程中，網(wǎng)格的數(shù)量決定了計(jì)算的精度，因此在進(jìn)行計(jì)算前，本文先驗(yàn)證了網(wǎng)格數(shù)量對(duì)計(jì)算精度的影響。分別選取了網(wǎng)格數(shù)量225671、446242、682542進(jìn)行計(jì)算精度驗(yàn)證，表3給出了3種網(wǎng)格數(shù)量下，計(jì)算穩(wěn)定時(shí)水的飽和度。計(jì)算結(jié)果表明，網(wǎng)格數(shù)量較小時(shí)，計(jì)算誤差較大，但網(wǎng)格數(shù)量提高大約3倍時(shí)，以網(wǎng)格數(shù)量最多的計(jì)算結(jié)果作為參考，中等網(wǎng)格數(shù)量相對(duì)誤差為4.81%。因此從計(jì)算精度和計(jì)算資源考慮，本次選中等網(wǎng)格數(shù)量，即446242。

表3 網(wǎng)格無關(guān)性驗(yàn)證

2.1.2 VOF模型可行性驗(yàn)證

為了驗(yàn)證VOF模型計(jì)算Sc-CO2的可行性。分別使用單相層流模型和VOF層流模型（Sc-CO2體積分?jǐn)?shù)為1）進(jìn)行驗(yàn)證，入口使用質(zhì)量流量入口，入口流量大小都為0.0005g/s，出口使用壓力出口[14]。兩種模型對(duì)應(yīng)的進(jìn)出口壓降相差為3.62×10–5Pa，由此驗(yàn)證VOF模型可用于計(jì)算Sc-CO2驅(qū)替含水多孔介質(zhì)。

2.2 Sc-CO2驅(qū)替

2.2.1 Sc-CO2驅(qū)替過程兩相界面分析

整個(gè)驅(qū)替過程中，Sc-CO2以恒定質(zhì)量流量從左側(cè)入口進(jìn)入，初始溫度設(shè)定為333K，環(huán)境壓力為15MPa，表面張力系數(shù)設(shè)定為0.045N/m，圖6給出了兩相界面與壁面接觸角分別為30°和60°時(shí)整個(gè)驅(qū)替過程中不同時(shí)刻的兩相分布圖。計(jì)算結(jié)果表明，在t=0.0005s時(shí)刻，Sc-CO2運(yùn)移到第一個(gè)對(duì)稱顆粒，兩相界面對(duì)稱分布；在t=0.001s，Sc-CO2運(yùn)移到第一個(gè)上下非對(duì)稱的小顆粒處，30°與60°接觸角出現(xiàn)了不同結(jié)果，30°接觸角選擇小孔隙吼道進(jìn)行突破，而60°接觸角選擇大孔隙吼道進(jìn)行突破。這是因?yàn)楫?dāng)吼道較小時(shí)，60°接觸角中吼道較小的地方形成毛細(xì)壓力，故對(duì)于非潤(rùn)濕相驅(qū)替潤(rùn)濕相，毛細(xì)壓力變成阻力，而相比于較小的接觸角時(shí)，毛細(xì)壓力成為驅(qū)替動(dòng)力，因此驅(qū)替過程中出現(xiàn)不同現(xiàn)象。隨著Sc-CO2連續(xù)的注入到多孔介質(zhì)中，在壁面形成的潤(rùn)濕相液膜也逐漸被驅(qū)替，在t=0.0027s時(shí)刻通過最后一個(gè)小顆粒。

2.2.2 毛細(xì)管數(shù)的影響

圖6 不同時(shí)刻驅(qū)水兩相分布圖

由于多孔介質(zhì)中含有不均勻的孔隙，因此為了探討孔隙大小對(duì)驅(qū)替過程的影響，本文通過改變Sc-CO2的注入流量來獲取不同的毛細(xì)管數(shù)。圖7給出了數(shù)值模擬結(jié)果與WANG等[7]在2012年通過實(shí)驗(yàn)獲取毛細(xì)管數(shù)影響Sc-CO2飽和度數(shù)據(jù)圖。計(jì)算結(jié)果表明，采用數(shù)值模擬下的不同毛細(xì)管數(shù)的Sc-CO2飽和度能很好的與實(shí)驗(yàn)數(shù)據(jù)吻合；隨著注入流量的增加，Sc-CO2飽和度下降，但隨著流量的再次加大，毛細(xì)管數(shù)達(dá)到Sc-CO2的飽和度變化不再明顯，最后趨于平穩(wěn)。

圖7 Sc-CO2飽和度與毛細(xì)管數(shù)關(guān)系

2.2.3 Sc-CO2注射壓力的影響

咸水層地質(zhì)封存深度大約為1000～2500m，不同深度下的封存過程中的環(huán)境壓力會(huì)有所不同。因此本文通過改變模擬中的環(huán)境壓力來實(shí)現(xiàn)不同封存壓力，壓力方案A、B、C分別選取10MPa、15MPa、20MPa。圖8給出了在t=0.0012s時(shí)刻Sc-CO2驅(qū)水過程兩相分布圖；3種方案殘余水飽和度隨驅(qū)替時(shí)間發(fā)生的變化關(guān)系如圖9所示。數(shù)值模擬結(jié)果表明，水的飽和度在10%以下時(shí)，3種方案所花費(fèi)的時(shí)間大致為0.0015s、0.0026s、0.0033s，因此封存壓力越低，驅(qū)替速率越快，其次由最后驅(qū)替達(dá)到穩(wěn)定狀態(tài)結(jié)果表明，環(huán)境壓力越低，殘余水飽和度越低。這是因?yàn)榈刭|(zhì)封存選取深度越高，即驅(qū)替環(huán)境壓力越高，這時(shí)Sc-CO2的密度越大，在同樣的質(zhì)量流量下，體積流量反而越小，因此Sc-CO2樣能運(yùn)移到更遠(yuǎn)處；同時(shí)驅(qū)替所遇到的阻力越大，導(dǎo)致驅(qū)替速率下降，但是驅(qū)替過程所花費(fèi)的時(shí)間并不與壓力變化成正比，壓力越高驅(qū)替速率變化減小。

圖8 t=0.0012s時(shí)刻兩相分布圖

圖9 H2O飽和度與驅(qū)替時(shí)間關(guān)系

2.2.4 Sc-CO2注射溫度的影響

由于地質(zhì)封存選取地層深度約為1000～2500m，該地層的環(huán)境溫度大致為313～383K，本文研究選取含水多孔介質(zhì)溫度為313K，采用不同注射溫度下的Sc-CO2來驅(qū)替，方案D、E、F溫度分別為313K、323K、333K。圖10給出了Sc-CO2在3種溫度條件、t=0.002s時(shí)刻的兩相分布圖；圖11給出了3種溫度條件下殘余水飽和度隨驅(qū)替過程的變化圖。計(jì)算結(jié)果表明，Sc-CO2所選取的溫度越高，驅(qū)替速率越快，在驅(qū)替中期，驅(qū)替速率大致和Sc-CO2溫度成線性關(guān)系，這是因?yàn)殡S著溫度增加，Sc-CO2黏度越小，流動(dòng)性強(qiáng)，即在相同時(shí)刻，方案F中Sc-CO2運(yùn)移到更遠(yuǎn)處，但是隨著時(shí)間的推移，驅(qū)替穩(wěn)定后，驅(qū)替效率大致相等。

2.2.5 表面張力的影響

圖10 t=0.002s時(shí)刻兩相分布圖

圖11 H2O飽和度與驅(qū)替時(shí)間關(guān)系

圖12 t=0.0018s時(shí)刻兩相分布圖

圖13 H2O飽和度與驅(qū)替時(shí)間關(guān)系

為了研究表面張力對(duì)兩相驅(qū)替過程的影響，方案G、I分別選取表面張力為0.03N/m、0.06N/m來進(jìn)行研究。圖12給出了兩種方案在t=0.0018s時(shí)刻兩相分布圖，圖13給出了殘余水飽和度隨驅(qū)替時(shí)間發(fā)生變化關(guān)系圖。計(jì)算結(jié)果表明，隨著表面張力的增加，驅(qū)替速率下降，同時(shí)驅(qū)替效率也下降。這是因?yàn)殡S著表面張力的減小，Sc-CO2在驅(qū)替過程中兩相毛細(xì)壓力越小，因此在相同的時(shí)間內(nèi)，圖12方案G中Sc-CO2能運(yùn)移到更遠(yuǎn)處；同時(shí)由于表面張力的增加，潤(rùn)濕相更易黏附在固體壁面，導(dǎo)致驅(qū)替效率有所下降，故殘余水飽和度更大。

2.2.6 接觸角的影響

不同的地質(zhì)條件下，巖石的浸潤(rùn)性會(huì)有所不同，也即表現(xiàn)為三相（巖石、水、Sc-CO2）接觸角不一樣，因此本文以此為研究出發(fā)點(diǎn)，選取不同的接觸角（方案J、K、M的接觸角分別為30°、60°、90°），在其他條件不變的情況下進(jìn)行數(shù)值模擬計(jì)算。3種不同接觸角方案在t=0.0016s時(shí)兩相分布見圖14；圖15給出了驅(qū)替過程中殘余水飽和度隨時(shí)間變化的過程。數(shù)值模擬結(jié)果表明，接觸角越小，驅(qū)替速率越快，但驅(qū)替效率有所下降。這是因?yàn)楫?dāng)巖石壁面表現(xiàn)為親水性時(shí)，隨著接觸角的減小，潤(rùn)濕相潤(rùn)濕性越好，這樣會(huì)在壁面形成液膜，使得非潤(rùn)濕相Sc-CO2在驅(qū)替過程中阻力減小，也即表現(xiàn)為方案J中Sc-CO2運(yùn)移到更遠(yuǎn)處；但是接觸角越小，潤(rùn)濕相在壁面形成液膜越厚，Sc-CO2驅(qū)替效率下降。

圖14 t=0.0016s時(shí)刻兩相分布圖

圖15 H2O飽和度與驅(qū)替時(shí)間關(guān)系

3 結(jié)論

本文建立了孔隙尺度下各向同性的三維多孔介質(zhì)模型，使用數(shù)值模擬方法研究了Sc-CO2注入含水多孔介質(zhì)的兩相驅(qū)替過程。分析了毛細(xì)管數(shù)、地質(zhì)封存壓力、Sc-CO2注入溫度、兩相表面張力、接觸角影響因素對(duì)Sc-CO2驅(qū)替含水多孔介質(zhì)中的運(yùn)移規(guī)律，得到如下結(jié)論。

毛細(xì)管數(shù)越大，Sc-CO2飽和度先下降，但隨著毛細(xì)管數(shù)對(duì)數(shù)進(jìn)一步加大到–5.5時(shí)，Sc-CO2的飽和度趨于平穩(wěn)；地質(zhì)封存地點(diǎn)越深，達(dá)到同一殘余水飽和度驅(qū)替所花費(fèi)的時(shí)間越長(zhǎng)，但花費(fèi)的時(shí)間并不與封存壓力成線性關(guān)系，當(dāng)壓力大于15MPa時(shí)，驅(qū)替速率大大下降；Sc-CO2注射溫度越高，驅(qū)替速率越快；隨著兩相表面張力系數(shù)的增加，驅(qū)替速率減小，同時(shí)驅(qū)替效率也越低；在壁面表現(xiàn)為潤(rùn)濕性條件下，接觸角越小，驅(qū)替速率越快，同時(shí)較小的接觸角，將導(dǎo)致過大的毛細(xì)壓力，這可能成為驅(qū)替的動(dòng)力。

[1] BODE S，JUNG M. Carbon dioxide capture and storage—liability for non-permanence under the UNFCCC[J]. International Environmental Agreements: Politics，Law and Economics，2006，6（2）：173-186.

[2] BENTHAM M，KIRBY M. CO2storage in saline aquifers[J]. Oil &Gas Science and Technology，2005，60（3）：559-567.

[3] CHADWICK R，ARTS R J，BERNSTONE C，et al. Best practice for the storage of CO2in saline aquifers-observations and guidelines form the SACS and CO2STORE projects[M]. British Geological Survey，2008.

[4] THOMAS S. Enhanced oil recovery- an overview[J]. Oil & Gas Science and Technology-Revue de l'IFP，2008，63（1）：9-19.

[5] BACHU S，SHAW J C，PEARSON R M. Estimation of oil recovery and CO2storage capacity in CO2EOR incorporating the effect of underlying aquifers[C]//Society of Petroleum Engineers，2004.

[6] SUEKANE T，SOUKAWA S，IWATANI S，et al. Behavior of supercritical CO2injected into porous media containing water[J].Energy，2005，30：2370–2382.

[7] WANG Y，ZHANG C Y，WEI N，et al. Experimental study of crossover from capillary to viscous fingering for supercritical CO2-water displacement in a homogeneous pore network[J].Environmental Science & Technology，2013，47：212-218.

[8] TREVISAN L，PINI R，CIHAN A，et al. Experimental investigation of supercritical CO2trapping mechanisms at the intermediate laboratory scale in well-defined heterogeneous porous media [J].Energy Procedia，2014，63：5646–5653.

[9] TREVISAN L，CIHAN A，F(xiàn)AGERLUND F，et al. Investigation of mechanisms of supercritical CO2trapping in deep saline reservoirs using surrogate fluids at ambient laboratory conditions[J].International Journal of Greenhouse Gas Control，2014，29：35–49.

[10] SASAKI K，F(xiàn)UJI T，NIIBORI Y，et al. Numerical simulation of supercritical CO2injection into subsurface rock masses[J]. Energy Conversion and Management，2008，49：54–61.

[11] CHOI H S，PARK H C，HUH C，et al. Numerical simulation of fluid flow and heat transfer of supercritical CO2in micro-porous media[J].Energy Procedia，2011，4：3786–3793.

[12] 馬瑾，胥蕊娜，羅庶，等. 超臨界壓力CO2在深部咸水層中運(yùn)移規(guī)律研究[J]. 工程熱物理學(xué)報(bào)，2012，33（11）：1971-1975.MA J，XU R，LUO S，et al. Core-scale experimental study on supercritical-pressure CO2migration mechanism during CO2geological storage in deep saline aquifers[J]. Journal of Engineering Thermophysics，2012，33（11）：1971-1975.

[13] XU R，LUO S，JIANG P X. Pore scale numerical simulation of supercritical CO2injection into porous media containing water[J].Energy Procedia，2011，4：4418-4424.

[14] 羅庶，胥蕊娜，姜培學(xué). CO2地質(zhì)封存突破壓力梯度的孔隙尺度數(shù)值模擬[J]. 工程熱物理學(xué)報(bào)，2011，32（5）：819-823.LUO S，XU R N，JIANG P X. Pore-sclae numerical simulation on breakthrough pressure gradient during CO2geological storage[J].Journal of Engineering Thermophysics，2011，32（5）：819-823.

[15] 蔣蘭蘭. CO2地質(zhì)封存多孔介質(zhì)內(nèi)氣液兩相滲流特性研究[D]. 大連：大連理工大學(xué)，2014.JIANG L L. Study on gas/liquid flow characteristic in porous media during CO2geological storage[D]. Dalian: Dalian University of Technology，2014.

[16] XU R N，JIANG P X. Numerical simulation of fluid flow in microporous media[J]. International Journal of Heat and Fluid Flow 2008，29（5）：1447-1455.

[17] POURYOUSSEFI S M，ZHANG Y W. Numerical investigation of chaotic flow in a 2D closed-loop pulsating heat pipe[J]. Applied Thermal Engineering，2016，98：617–627.

[18] CARCIOF B A M，PRAT M，LAURINDO J B. Homogeneous volume-of-fluid（VOF）model for simulating the imbibition in porous media saturated by gas[J]. Energy & Fuels，2011，25（5）：2267-2273.

[19] HIRT C W，NICHOLS B D. Volume of fluid （VOF） method for the dynamics of free boundaries[J]. Journal of Computational Physics，1981，39：201-225.

Pore-scale two-phase numerical simulation of supercritical carbon dioxide displacement of water

JIANG Shuisheng，ZHAO Wandong，ZHANG Ying，LI Peisheng，WANG Zhaotai，ZHONG Yuan
（School of Mechanical and Electrical Engineering，Nanchang University，Nanchang 330031，Jiangxi，China）

The computational method of pore-scale porous media is established by using Volume of Fluid（VOF）numerical simulation method based on the background of CO2storage in Saline aquifer，and the migration mechanism of supercritical carbon dioxide（Sc-CO2）injection into porous media containing water is studied. The effects of capillary number，geological storage pressure，Sc-CO2injection temperature，two-phase surface tension coefficient and contact angle on the two-phase migration rate and the displacement efficiency were analyzed. At the same time，displacement efficiency was compared with experimental data under different capillary numbers. The result showed that with the increase of capillary number，displacement efficiency decreases first and then becomes stable，and the numerical simulation of displacement efficiency agrees well with experimental data under different capillary numbers. Under the same porosity，when wall surface is hydrophilic，wall surface wettability is better. The faster the rate of displacement，while the displacement efficiency decreased. At the same time，the lower geological storage pressure，the higher the injection temperature.Small surface tension coefficient has a higher displacement rate and efficiency.

supercritical carbon dioxide；saline aquifer storage；pore-scale；numerical simulation；two-phase flow；capillary number；contact angle

TK124；TE122

A

1000–6613（2017）11–3955–08

10.16085/j.issn.1000-6613.2017-0567

2017-03-31；修改稿日期2017-06-22。

國(guó)家自然科學(xué)基金（51566012，11562011）、江西省科技廳支撐項(xiàng)目（2009BGA01800）及江西省研究生創(chuàng)新基金（YC2017－S056）項(xiàng)目。

姜水生（1957—），男，碩士，教授。聯(lián)系人張瑩，博士，教授，研究?jī)?nèi)容為復(fù)雜熱流場(chǎng)模擬。E-mail：yzhan2033@163.com。