數(shù)值模擬在砂巖型鈾礦流體運移機制研究中的應(yīng)用

韋曉艷，李增華，杜鵬飛，嚴兆彬，于炳飛，李勝苗，劉傳東

（1.東華理工大學 江西省放射性地學大數(shù)據(jù)技術(shù)工程實驗室，江西 南昌 330013；2.東華理工大學 地球科學學院，江西 南昌 330013；3.核工業(yè)二〇八大隊，內(nèi)蒙古 包頭 014010；4.中國冶金地質(zhì)總局中南地質(zhì)調(diào)查院，湖北 武漢 430080；5.湖南省地質(zhì)調(diào)查院，湖南 長沙 410016）

砂巖型鈾礦是指產(chǎn)于近地表砂巖中（部分為礫巖、粉砂巖、泥巖）以外生成礦作用為主形成的鈾礦床［1］。近年來，隨著地浸開采技術(shù)的不斷完善，砂巖型鈾礦具有開采成本低且相對環(huán)保等優(yōu)點，是各個國家現(xiàn)今重點開發(fā)對象和研究熱點［2］。砂巖型鈾礦的產(chǎn)出都與沉積盆地有關(guān)，其形成為多種地質(zhì)因素綜合控制的結(jié)果，其中流體在成礦作用過程中扮演著十分重要的角色，是鈾活化、遷移和富集的主要介質(zhì)［3-4］。流體運移機制的研究是建立綜合成礦模式的重要組成部分［5］。

前人研究表明，砂巖型鈾礦的成礦流體運移驅(qū)動力主要有4 種：地勢差、流體超壓、構(gòu)造變形和流體密度差［6-11］。流體驅(qū)動力既有可能以單一驅(qū)動力為主、也可能共同作用于成礦流體。在成礦作用過程中，若存在多源驅(qū)動力時，以哪種驅(qū)動力為主，以及它們之間如何平衡？目前，流體動力學數(shù)值模擬是解決這些問題最直接有效的技術(shù)手段，它可以重現(xiàn)并定量化的研究成礦流體的運移過程，以及界定有利于成礦的地質(zhì)條件［5］。國際上，流體動力學數(shù)值模擬在砂巖型鈾礦中已有廣泛的應(yīng)用［6-7］。我國在這方面的應(yīng)用研究還相對較少，僅Xue等（2011）［8］、郝偉林（2006）［9］等做了地勢差和流體超壓驅(qū)動方面的研究工作，研究基礎(chǔ)十分薄弱，亟需開展這方面的工作。本文在分析盆地流體運移的理論模型、控制因素的基礎(chǔ)上，總結(jié)砂巖型鈾礦中流體驅(qū)動機制的研究進展，通過實例分析來揭示流體動力學數(shù)值模擬在該領(lǐng)域的應(yīng)用前景，對今后研究方向提出展望，以資同行討論和參考。

1 盆地流體運移機制

1.1 盆地流體運移的基本理論

盆地流體是指存在于盆地沉積物孔隙空間并通過其進行流動的任何流體（包括烴類和水）［11］。盆地流體廣泛參與了沉積物的成巖、后生、成油和成礦過程。流體的流動通常遵循達西定律，可表示為：

式中，q—流速/（m·s-1）；k—滲透率/m2；μ—流體粘度/（kg·m-1·s-1）；L—滲流路徑長度/m；Φ—流體勢/Pa，其定義為：

式中，ρ—流體密度/（kg·m-3）；g—重力/（m·s-2）；z—高程/m；p—流體壓力/Pa。流體勢（Φ）是指地層中某一點的壓能與該點相對某基準面（如海平面）的位能之和。流體總是從高勢區(qū)向低勢區(qū)運移。

由公式（1）和（2）可以看出，流體的流動速度與巖石的滲透率（k）、流體的粘度（μ）、密度（ρ）、壓力（p）和高差（z）有關(guān)。

從地質(zhì)學的角度看，所有這些參數(shù)，尤其是流體壓力和密度、巖石形變和地形起伏，都與構(gòu)造環(huán)境密切相關(guān)。不同構(gòu)造環(huán)境具有不同的流體壓力背景、熱力場和地勢差，因而造成流體流動的驅(qū)動力不同。具體來看，包括：1）與地形起伏有關(guān)的地勢差驅(qū)動（又稱重力驅(qū)動）；2）流體超壓驅(qū)動；3）密度差驅(qū)動（又稱浮力驅(qū)動）；4）構(gòu)造變形。

1.2 地勢差（重力）驅(qū)動的盆地流體運移

地形驅(qū)動的流體流動亦稱重力驅(qū)動流動，是流體在靜水壓力體系下，由于地勢差的存在而驅(qū)動的流體流動。

在地形起伏區(qū)，大氣降水沿盆地邊緣出露的透水層向下滲入盆地，在重力驅(qū)動下，盆地流體呈連續(xù)流動并向較低部位排出。大規(guī)模流體的流動方向是從盆地邊緣流向盆地中心，從盆地的淺部到盆地的深部，又稱為“向心流”或“滲入型流體”［11-12］。

1.3 流體超壓驅(qū)動的盆地流體運移

流體超壓指的是在地下某一深度流體壓力大于靜水壓力的值［13］。流體超壓的產(chǎn)生是由于孔隙流體支撐部分或全部上覆巖層的壓力而不僅僅是上覆孔隙水的壓力（靜水壓力）。造成這一現(xiàn)象的原因與孔隙度的減少或流體體積的增加有關(guān)，當孔隙度變小或流體體積變大的速率大于流體流動的補償值時就會出現(xiàn)流體超壓現(xiàn)象。流體總是由高超壓區(qū)向低超壓區(qū)流動［5］。

在沉積盆地中，孔隙度的減小主要是由于不均衡沉積壓實作用。上覆沉積層的加厚而引起的壓實作用驅(qū)動的流體運移，是促進沉積物中孔隙水排出的主要因素。壓實作用下排出的水將向流體超壓較小的部位流動，即主要向上、向盆地邊緣、向盆地內(nèi)部的水下高地部位流動［10］。

盆地內(nèi)的生烴效應(yīng)也會引起流體超壓［10-11］。烴源巖在向油氣轉(zhuǎn)化過程中，固態(tài)有機質(zhì)會形成液態(tài)和氣態(tài)烴類和水，使其密度減小，體積膨大，烴源巖內(nèi)部的壓力因生成的流體進入到烴源巖的流體系統(tǒng)中而增大，與周圍巖層形成壓差，驅(qū)動流體從烴源巖中排出。流體體積增加的原因還可以是流體受熱膨脹或礦物脫水作用，結(jié)果同樣形成流體超壓驅(qū)動流體流動［5］。

1.4 密度差（浮力）驅(qū)動的盆地流體運移

公式（2）表明流體的密度變化也會影響流體勢的變化，因此也可以控制流體的流動。在沉積盆地中，流體的密度主要與溫度和鹽度有關(guān)。一般密度與溫度成反比，與鹽度成正比［5］。

溫度對流體密度的影響，主要有地溫梯度（溫度隨深度升高）和局部熱異常（如巖漿侵入體）兩種情況。這種由溫度變化導(dǎo)致的密度差而引起的流動方式為熱對流。理論條件下的自由熱對流主要由瑞利（Rayleigh）數(shù)表現(xiàn)，與地溫梯度、巖石滲透率、地層厚度等其他熱力學參數(shù)有關(guān)。

1.5 巖石形變驅(qū)動的盆地流體運移

巖石的形變以及產(chǎn)生的構(gòu)造一直被認為是地殼中流體流動的主要動力機制之一，巖石形變會造成局部或區(qū)域性的流體勢梯度，從而推動流體的流動。當巖石產(chǎn)生裂縫時，巖石的應(yīng)力得到釋放，巖石產(chǎn)生擴容，在裂縫區(qū)會產(chǎn)生一個負壓差從而驅(qū)動流體運移；當巖石上覆地層遭受剝蝕時，下伏地層巖石應(yīng)力得到釋放，巖石發(fā)生擴容產(chǎn)生負壓差從而驅(qū)動流體運移。

2 砂巖型鈾礦流體運移機制研究現(xiàn)狀

2.1 地勢差驅(qū)動為主導(dǎo)的流體運移

眾多國內(nèi)外研究實例表明，地勢差是砂巖型鈾礦床盆地流體運移的重要機制之一。如Sanford（1990）［6］認為美國San Juan 盆地板狀砂巖型鈾礦床的形成與地勢差驅(qū)動的流體流動有關(guān)，其中盆地內(nèi)的Grants 鈾礦田主要產(chǎn)于流體向上流動的排水區(qū)。Sanford（1994）［7］進一步利用MODFLOW 軟件對San Juan 盆地的地下水流系統(tǒng)進行了模擬，結(jié)果顯示進入地下水系統(tǒng)的大氣水主要受重力驅(qū)動流向盆地東部和北部。此外，地勢差驅(qū)動的流體運移距離可達幾百千米到上千千米。如Banner 等（1989）［14］對美國密蘇里中部泉水和自流井含鹽地下水的同位素和微量元素研究表明，這些流體最可能起源于西部約1 000 km 科羅拉多落磯山山前帶的大氣降水補給。Bentley 等（1986）［15］在澳大利亞大自流盆地（Great Artesian Basin）的研究中，利用同位素示蹤方法也揭示出大規(guī)模的深部重力驅(qū)動流體系統(tǒng)。

國內(nèi)學者也有類似的認識。李啟榮（2000）［16］對新疆布爾津盆地砂巖型鈾礦成礦水文地質(zhì)條件進行分析研究，提出地勢差對流體運移具有驅(qū)動作用，由高水位的補給區(qū)向低水位的排泄區(qū)運動。許來生（2005）［17］對內(nèi)蒙古武川-烏蘭花盆地砂巖型鈾礦區(qū)進行了水文地質(zhì)條件分析，識別出不同水文地質(zhì)單元內(nèi)的重力驅(qū)動流體運移。郝偉林（2006）［9］在水文地質(zhì)和同位素水文學研究的基礎(chǔ)上，用MODFLOW 研究了鄂爾多斯盆地東勝地區(qū)的古水動力系統(tǒng)，結(jié)果表明地勢的高差控制了流體的運移。

Xue 等（2010）［18］利用Basin2 軟件對鄂爾多斯盆地內(nèi)流體流動系統(tǒng)進行數(shù)值模擬并建立了成礦流體運移模型，結(jié)果顯示盆地流體流動主要受地勢差控制，在盆地邊緣相對抬升和盆地內(nèi)地層輕微變傾過程中，地勢差驅(qū)動流體從盆地邊緣向下流向盆地中心。廉康等（2020）［19］對柴達木盆地北緣魚卡地區(qū)砂巖型鈾成礦條件分析指出魚卡地區(qū)地形特征北高南低、東西兩端高，中部低，區(qū)內(nèi)淺層水受地勢差驅(qū)動由北山區(qū)向盆地中側(cè)向補給。

2.2 沉積壓實和生烴作用導(dǎo)致的流體超壓

流體超壓可能是砂巖型鈾礦床流體運移的又一重要機制。許多研究表明，不均衡沉積壓實作用和含油盆地內(nèi)的生烴作用會造成流體超壓的現(xiàn)象。其中沉積壓實作用驅(qū)動沉積物間的封存水、孔隙水或裂隙水向盆地壓力較小的覆蓋層較薄區(qū)、盆地邊緣區(qū)、斷裂發(fā)育區(qū)和斷隆區(qū)流動泄壓排泄，對于鈾礦的形成具有一定的促進作用［10，12］。

但沉積壓實作用引起的流體流動跟其他驅(qū)動力相比一般較弱，不占主導(dǎo)地位。如Sanford（1994）［7］通過模擬San Juan 盆地的地下水流系統(tǒng)，顯示進入地下水系統(tǒng)的大氣水主要受重力驅(qū)動，沉積壓實起到的作用很小。Bj?rlykke（1994）［20］針 對沉積 盆地流體流動系統(tǒng)進行模擬，結(jié)果顯示與重力驅(qū)動的水流相比，壓實水是地下水流量的一個非常小的組成部分，由壓實驅(qū)動的孔隙水流速通常比重力驅(qū)動流速度低幾個數(shù)量級。Xue 等（2010）［18］利用Basin2 對鄂爾多斯盆地流體流動進行了模擬，模擬過程中只考慮不均衡沉積壓實作用時，盆地上覆沉積壓實作用產(chǎn)生的超壓力驅(qū)動滲出流，從盆地中心和下部流向盆地邊緣和淺部。當同時考慮沉積壓實與地勢差時，沉積壓實驅(qū)動的流體流動速率較為緩慢，因而流體流動主要還是受地勢差控制。但是，在沉積壓實作用的基礎(chǔ)上考慮盆地生烴效應(yīng)對流體流動的影響時，流體超壓比沒有考慮生烴效應(yīng)情況下的流體超壓值大了許多，并驅(qū)動流體向上、向邊部流動。生烴效應(yīng)是流體超壓產(chǎn)生的主要貢獻［18，21］。

2.3 構(gòu)造對流體運移的作用

盆地的構(gòu)造背景及其演化對于大規(guī)模的成礦作用具有重要的控制作用。盆地沉積后在成礦階段，構(gòu)造直接影響著成礦流體的補-徑-排體系，制約著鈾礦的產(chǎn)出［22］。一方面構(gòu)造是流體運移的重要通道，同時構(gòu)造形變也是驅(qū)動盆地流體運移的重要機制。

李啟榮（2000）［16］分析了新疆布爾津盆地砂巖型鈾礦的水文地質(zhì)條件，認為滲出承壓系統(tǒng)水動力區(qū)（帶）沿盆內(nèi)隱伏深大斷裂兩側(cè)成狹窄條帶狀分布，斷層脈狀承壓水在構(gòu)造擠壓應(yīng)力作用下沿斷裂自地殼深部往上運移至沉積蓋層，沿層間集水砂體向斷裂兩側(cè)運動，并與地勢差驅(qū)動的層間滲入水相遇，二者的交替混合作用導(dǎo)致鈾的沉淀并富集于層狀砂體中。

另外深部熱流體對成礦的作用得到研究者們的關(guān)注［23-24］。林雙幸等（2017）［12］認為哈薩克斯坦、烏茲別克斯坦、俄羅斯（外烏拉圭地區(qū)）以及我國北方的中新生代盆地的鈾礦床均受構(gòu)造、巖漿和火山活化作用影響，受淺成中基性巖漿活動和貫穿性斷裂作用，攜帶成礦元素的深部熱流體沿活動斷裂和羽狀裂隙貫入到透水層并向排泄區(qū)漫流，發(fā)生熱流體蝕變等地化作用形成鈾礦。

開魯盆地錢家店-白興吐礦床及其周圍斷裂構(gòu)造十分發(fā)育，尤其是著名的西拉木倫河斷裂在礦床中直穿而過，為流體的運移提供了通道。李建國等（2020）［25］通過對錢家店砂巖型鈾礦床進行蝕變礦物填圖，認為在晚白堊世末期持續(xù)隆升的構(gòu)造背景下，沿斷裂上升的酸性還原性流體對成礦貢獻巨大。

近年的鈾礦年代學研究也表明不同時期的盆地構(gòu)造隆升、掀斜等運動與成礦具有時空上的密切聯(lián)系［26-27］，這些構(gòu)造運動可能為含氧含鈾流體進入成礦目的層提供了動力。

2.4 流體密度差（浮力）驅(qū)動流體運移

流體的密度差同樣能夠驅(qū)動流體的運移，影響流體的流向以及排泄區(qū)和補給區(qū)。Zhang等（2018）［26］對鄂爾多斯盆地北部東勝鈾礦的研究認為，流入直羅組的地下水以大氣降水為主，其密度遠小于沉積物中已經(jīng)存在的鹵水的密度，這種密度差會導(dǎo)致大氣水集中在層序上部，進而控制著鈾礦化的位置。此外，Sanford（1994）［7］利用MODFLOW 軟件對美國San Juan盆地板狀砂巖型鈾礦床的地下水流系統(tǒng)進行模擬時研究了鹽湖水和淡水密度差對流體流動的影響，發(fā)現(xiàn)在向下流動的淡水流體和從鹽湖向外流動的堿性流體之間存在一個界面，界面近水平并向下延伸，但實際的流動方向及速率取決于鹽湖水密度和地下孔隙水的密度等因素的變化。

3 實例分析

作為砂巖型鈾礦的成礦流體運移機制研究的一次初步嘗試，本次實例模擬采用FLAC3D 軟件研究砂巖盆地受擠壓后在抬升過程中壓力的變化以及流體流動情況。

3.1 使用軟件簡介

本次實例模擬采用FLAC 軟件（Fast Lagrangian Analysis of Continua，包 括2D 和3D），該軟件是美國Itasca 公司開發(fā)的，可較好地模擬巖石的構(gòu)造變形，還可進行與熱力場和流體場的耦合研究，在礦床學中常用來研究構(gòu)造變形與流體運移［28］。FLAC 軟件采用的是顯式的有限差分法，以節(jié)點為計算對象，將力和質(zhì)量均集中在節(jié)點上，然后通過運動方程在時域內(nèi)進行求解，對于每個節(jié)點，由應(yīng)力及外力利用虛功原理求節(jié)點不平衡力，再由不平衡力求出節(jié)點速率；對于每個單元，由節(jié)點速率求出應(yīng)變增量，再由應(yīng)變增量求應(yīng)力增量及總應(yīng)力。

FLAC 軟件適用于非線性、大位移、大應(yīng)變問題，無需額外復(fù)雜的計算。FLAC 的工作流程為建立基礎(chǔ)的幾何模型，定義模型類型（本構(gòu)關(guān)系），輸入各種地層參數(shù)（如粘聚力，抗拉強度和摩擦力等），確定初始條件以及邊界條件，然后對工作區(qū)的構(gòu)造應(yīng)力場、熱力場以及流體場進行模擬，給出模擬結(jié)果。

3.2 建立地質(zhì)模型

針對“層間氧化帶成礦模式”，本次數(shù)值模擬整體采用的是泥-砂-泥結(jié)構(gòu)，構(gòu)建的幾何模型寬10 km，高5 km，模型包含盆地蓋層和基底兩大部分。蓋層從下至上分別為：泥巖層→粉砂巖層→砂質(zhì)泥巖層→砂巖層→泥巖層（圖1）。

圖1 用于數(shù)值模擬的地質(zhì)概念模型Fig.1 Conceptual model used for numerical simulation

3.3 巖石物理參數(shù)

模擬采用的是FLAC3D 內(nèi)置的摩爾-庫倫本構(gòu)模型。各地質(zhì)單元的滲透率和孔隙度（表1），主要根據(jù)鄂爾多斯盆地不同地質(zhì)單元的已有實測參數(shù)［29-30］，而力學參數(shù)主要來自于FLAC3D 手冊或相似巖性的測試數(shù)據(jù)［31-32］。

表1 數(shù)值模擬中巖石力學和流體參數(shù)Table 1 Parameters of different rock units used in the numerical simulation

3.4 初始條件、邊界條件

巖石的孔隙壓力開始狀態(tài)為靜水壓力，模型地表為透水邊界，其他邊界設(shè)為不透水邊界，隨后進行初始地應(yīng)力平衡得到平衡狀態(tài)。在模型左側(cè)1 km 的范圍內(nèi)的格子初始了位移速度（10-8m/s）來模擬擠壓和抬升的過程。

3.5 結(jié)果分析與討論

經(jīng)過4%的擠壓后（即左側(cè)邊界向右位移400 m），模型中各地層已發(fā)生明顯變形，出現(xiàn)明顯的隆起（圖2），隆起最大處達到420 m，地層呈現(xiàn)掀斜狀態(tài)（圖3）。在深度2 km 處的砂質(zhì)泥巖中出現(xiàn)流體超壓（圖4）。流體超壓的存在導(dǎo)致砂質(zhì)泥巖中的流體往下部流動，而地表由于形成了地勢差，驅(qū)動流體砂巖中的流體從左往右流動（圖5）。

圖2 模型擠壓4%后地層的變形情況Fig.2 Deformation of strata after 4% bulk shortening

圖3 模型擠壓4%后模型抬升位移分布Fig.3 Displacement of Z direction after 4% bulk shortening

圖4 模型擠壓4%后流體壓力分布，在2 km 深處的砂質(zhì)泥巖中有流體超壓出現(xiàn)Fig.4 Distribution of pore pressure after 4% bulk shortening and the overpressure within the sandy mudstone at depth of 2 km

圖5 模型擠壓4%后流體運移方向和速度Fig.5 Fluid direction and velocity after 4% bulk shortening

該擠壓模型揭示了盆地在抬升隆起過程中，對流體壓力和流動方向具有重要的控制作用。在盆地隆起過程中，成礦流體可能由于局部超壓的影響，流體容易沿滲透率大的砂巖中匯集，有利于礦化的發(fā)生。同時地表的隆起也造成地勢差，驅(qū)動流體的遷移。

4 結(jié)語與展望

砂巖型鈾礦的形成經(jīng)過了一個長期的、動態(tài)的、復(fù)雜的過程，既包括物理（動力）過程也包括化學（反應(yīng)）過程，流體在成礦過程中起到關(guān)鍵作用。因此，流體動力學的研究不僅對建立成礦模型至關(guān)重要，對于礦床的空間定位和礦體的控制因素也有一定的指導(dǎo)意義。目前對于流體性質(zhì)、來源的研究較為普遍，但對于不同來源流體的驅(qū)動力、壓力狀態(tài)、流速、流動方向和規(guī)模等物理過程的研究還較為缺乏。成礦流體動力學數(shù)值模擬技術(shù)是解決這一問題的重要途徑之一，其研究結(jié)果可以提高對礦床成礦過程的認識，能合理地再現(xiàn)復(fù)雜的流體運移過程，定量地分析不同地質(zhì)因素對流體運移的影響。通過與已有地質(zhì)研究成果進行相互查證和綜合對比，模擬結(jié)果可用來較為準確地揭示成礦流體的運移規(guī)律，并可直接或間接地應(yīng)用于找礦勘探。因此在砂巖型鈾礦成礦規(guī)律和成礦預(yù)測方面，應(yīng)加強流體動力學的研究。