潘保芝,張 瑞,劉 坤,李 丁,蔣必辭

(1.吉林大學地球探測科學與技術學院,吉林長春130026；2.中國石油大學(華東)地球科學與技術學院,山東青島266580)

為取代傳統密度測井中的同位素中子源,各大石油公司都在研究利用脈沖中子源測量地層密度的方法。脈沖中子源是一個小型的加速器,將氘核加速轟擊鋯氚靶上的氚核,通過聚變反應產生快中子[1]。脈沖中子源以一定的脈沖寬度和重復周期向地層發(fā)射能量為14MeV的中子束。與常規(guī)中子源相比,脈沖中子源具有源強度高且質量輕、釋放出的中子能量單一等特點,特別是能夠依賴電信號以脈沖的形式向外界發(fā)射中子,而當電信號停止,中子便不再發(fā)射,故沒有放射性危害。

關于脈沖中子源密度測量的研究,已經取得一定的成就。1997年,康普樂公司首先研制出了利用快中子非彈散射伽馬射線獲取地層密度的測井儀器[1]；2005年,斯倫貝謝公司推出了包括脈沖中子密度測井在內的隨鉆測井平臺[2-3]；2013年Atfeh等[4]利用脈沖中子源密度測井方法對大規(guī)模碳酸鹽巖儲層進行了實例研究,進一步證明了脈沖中子源方法進行地層評價的可用性。

前人對于脈沖中子源密度測井的研究,大多以非彈散射伽馬射線為研究對象,目前相關理論已基本成熟。本文的研究對象是俘獲伽馬射線,將俘獲伽馬射線作為密度測井的伽馬源進行研究,模擬分析經地層吸收后剩余伽馬射線強度與地層密度以及探測源距的關系,從而指導脈沖中子源密度測井儀器的設計。

蒙特卡洛數值模擬方法是模擬微觀粒子運輸的重要工具,在核測井方面有很大的應用價值。它成本低,周期短,并在一定程度上彌補了無法獲得實體物理模型的缺陷[5]。我們采用蒙特卡洛數值模擬方法模擬并記錄粒子的軌跡,得到中子與原子核碰撞以及擴散后的位置,研究中子被俘獲產生的俘獲伽馬射線強度的分布,分析經地層吸收后剩余伽馬射線強度與地層密度的關系,并討論測量條件對剩余伽馬射線強度分布的影響。

1 模擬原理及建立模型

脈沖中子源發(fā)射的快中子進入地層后,依次發(fā)生減速作用損失能量以及擴散作用后被俘獲[6]。發(fā)生減速作用和被俘獲時分別產生非彈散射伽馬射線和俘獲伽馬射線。我們將產生的俘獲伽馬射線作為密度測井的伽馬源,該伽馬射線進入地層后,會被地層衰減,部分被吸收,通過模擬經地層吸收后剩余伽馬射線的強度,便可以得出其強度與地層密度的關系。

首先建立一個各向同性均勻二維介質平面模型,為了使效果更直觀清楚,采用直角坐標系研究極坐標問題,將中子源設定在坐標原點(0,0),模擬范圍為200cm×200cm,骨架是石英,孔隙度(φ)可變,孔隙流體是含氯(Cl)的礦化水。脈沖中子源發(fā)射大量的快中子,假設快中子進入地層時的初始能量為14MeV,不考慮最初的非彈散射,快中子與孔隙中的流體或者巖石的骨架發(fā)生彈性碰撞,其能量逐漸減弱,變成熱中子,之后經過擴散作用,最后被地層中的核素俘獲并產生俘獲伽馬射線。大量的俘獲伽馬射線進入地層之后,經康-吳效應和光電效應將被地層吸收一部分,最終得到剩余伽馬射線的分布。圖1為用俘獲伽馬射線的衰減探測地層密度示意圖。

圖1 用俘獲伽馬射線的衰減探測地層密度示意圖解

2 中子與物質作用過程的模擬

2.1 減速過程的模擬

脈沖中子源發(fā)射的快中子,依次經過地層的減速作用和擴散作用而被地層核素俘獲。首先模擬地層對快中子的減速作用過程。不考慮最開始的非彈散射,只模擬彈性散射減速作用。而減速作用主要是地層中的氫核素(H)引起的,所以只考慮與孔隙中H的彈性碰撞過程。快中子彈性碰撞過程中能量的變化為

(1)

式中：E0表示碰撞前的能量；E′表示碰撞后的能量；M為靶核的質量；m為快中子的質量；θ為彈性散射角度。

對于大量快中子的模擬,建立分別表示x坐標、y坐標和能量E的3個二維數組xij,yij和Eij,數組下標表示第i個快中子的第j次碰撞；為了獲得足夠數量的計算,用new(·)函數將數組定義為動態(tài)的,用rand(·)函數產生隨機數賦給每一個散射角度,兩次碰撞之間的路程為1cm[7],初始能量為14MeV,終止能量為熱中子能量0.025eV。

模擬了中子數量N=1000時不同孔隙度情況下熱中子的分布(圖2),圖2a,圖2b和圖2c分別是地層孔隙度為5%,40%,100%(水)時1000個熱中子的初始位置云圖。從圖2中可以看出,熱中子是圍繞中子源呈圓形分布的,在離源近的地方分布較密集,離源遠的地方分布較稀疏；熱中子云的半徑隨介質孔隙度的增大而逐漸減小,這是因為在減速過程中,H是最重要核素,含氫指數越大,快中子越容易被減速。由于H主要出現在油或水中,所以在地層中不含泥或不含結晶水的情況下,可以通過脈沖中子測井研究地層中的含氫量,從而求取孔隙度。通過熱中子云的半徑來判斷介質的減速長度,由圖2c可得水的減速長度約為7.5cm,這剛好與實際查得的資料相符,實際查閱水的減速長度為7.7cm[6]。

圖2 不同地層孔隙度條件下熱中子的初始位置云圖(N=1000)

2.2 擴散過程和俘獲過程的模擬

(2)

式中：∑a與∑s分別為介質的宏觀俘獲截面和宏觀散射截面。已知水的∑a=0.0220cm-1，∑s=2.6800cm-1；石英的∑a=0.0034cm-1，∑s=0.2680cm-1[6]。在得到每個熱中子的初始位置之后,熱中子要在介質中擴散,擴散的距離由(2)式計算。擴散的角度隨機,用rand(·)函數生成1～100的隨機數,當隨機數小于孔隙度時,認為熱中子在水中擴散,當隨機數大于孔隙度時,認為熱中子在石英中擴散。圖3和圖4分別顯示了孔隙度為20%和50%時擴散前、后的1000個熱中子位置云圖。

由圖3和圖4可以看出,熱中子擴散使得中子云半徑增大,孔隙度越大,熱中子擴散長度越小。這是因為當孔隙度增大時,巖石骨架含量相應減少,骨架的熱中子擴散長度比孔隙水的熱中子擴散長度大,所以綜合結果導致熱中子擴散長度減小。

圖3 地層孔隙度為20%時擴散前(a)、后(b)的熱中子位置云圖(N=1000)

圖4 地層孔隙度為50%時擴散前(a)、后(b)的熱中子位置云圖(N=1000)

熱中子擴散結束即被核素俘獲同時產生俘獲伽馬射線,對于幾種常見的核素來說,除了H核素之外,氯(Cl)、硅(Si)、鋁(Al)、鈣(Ca)俘獲熱中子產生的伽馬射線能量相差并不大,都在7MeV左右[6]。所以模擬時將放出的伽馬射線能量定為7MeV。這些俘獲伽馬射線便可作為密度測井的“伽馬源”。

在這些核素中,氧(O)的俘獲截面與其它核素相比低很多[6],故將O的俘獲作用省略,只模擬H,Cl和Si。依然模擬一個巖石骨架為石英的純地層,而孔隙流體為具有一定礦化度的水,將被3種不同核素俘獲的熱中子位置分別進行統計,得到其熱中子俘獲位置云圖。圖5是孔隙度為20%,氯含量占孔隙體積的百分比(VCl)為10%時,1000個熱中子被Cl,Si,H俘獲的位置云圖。圖5中被Cl俘獲放出的伽馬射線最多,Si次之,而H最少,這是因為三者的俘獲截面不同以及含量不同所引起的。

圖5 熱中子被Cl,Si,H俘獲的位置云圖(φ=20%,%,N=1000)

3 俘獲伽馬射線吸收過程的模擬

3.1 俘獲伽馬射線強度分布的模擬

已知熱中子被Cl和Si俘獲產生兩種能量的伽馬量子,被H俘獲產生一種伽馬量子[6],所以可以計算這個“伽馬源”的強度。將熱中子被俘獲的區(qū)域進行網格劃分,每個網格大小為5cm×5cm,求出每個單位面積內的伽馬量子數,統計得出俘獲伽馬射線的強度分布。為了找出更準確的統計規(guī)律,將中子數量增加到10000個進行模擬。圖6分別是孔隙度為5%,10%和20%時10000個熱中子產生的俘獲伽馬射線強度(Q)分布圖,這些也就是分散式的伽馬源。

由圖6可以看出,俘獲伽馬射線強度在中子源附近較大,隨著與中子源的距離L的增大,伽馬射線強度迅速減弱。并且,孔隙度小的地層伽馬射線強度分布較平緩,峰值較??；孔隙度大的地層伽馬射線強度分布較陡峭,峰值較大。即孔隙度越大,俘獲伽馬射線強度隨源距下降得越快。

圖6 不同孔隙度條件下10000個熱中子俘獲伽馬射線強度分布圖解

對圖6選取y=0,x>0的剖面進行繪圖,結果如圖7所示。從圖7中可以看出,俘獲伽馬射線強度隨源距L的增大而降低。當L>100cm時,俘獲伽馬射線強度已經趨近于零；當L≈35cm時,對于不同孔隙度的地層俘獲伽馬射線強度幾乎交于一點,此時俘獲伽馬射線強度與地層孔隙度無關；當L<35cm時,孔隙度高的地層俘獲伽馬射線強度大；當L>35cm時,孔隙度低的地層俘獲伽馬射線強度大。

圖7 10000個熱中子俘獲伽馬射線強度隨源距的變化

3.2 俘獲伽馬射線吸收過程的模擬

俘獲伽馬射線進入地層后,會被地層中的元素吸收,已知剩余伽馬射線強度J滿足如下關系式[6]：

(3)

式中：i表示面積單元序號;Qi表示第i個單元內俘獲放出的伽馬射線強度；ri表示第i個單元到探測點的距離；n表示單元的總個數；μm為質量吸收系數；ρ為地層密度。

模擬時首先進行網格劃分,共劃分成1600個5cm×5cm的網格,即n=1600。統計每個網格內的俘獲伽馬射線個數Qi,將質量吸收系數μm設定為0.17cm2/g[6]。假設地層孔隙度分別為5%,10%,20%,30%,40%,對于每一個孔隙度都分別假設地層密度為2.1,2.2,2.3,2.4,2.5g/cm3,統計在0,10,20,30,……,140cm源距處的剩余伽馬射線強度,從而得出其在縱向上的分布。圖8為剩余伽馬射線強度隨孔隙度變化(ρ=2.4g/cm3)和隨密度變化(φ=5%)的分布。

從圖8a可以看出,當L≈55cm時,剩余伽馬射線強度不受孔隙度的影響；當L<55cm時,剩余伽馬射線強度隨孔隙度的增大而增大；當L>55cm時,剩余伽馬射線強度隨孔隙度的增大而減小。從圖8b中可以看出,剩余伽馬射線強度隨密度的增高而減小。

圖8 10000個熱中子產生的剩余伽馬射線強度隨孔隙度變化(ρ=2.4g/cm3)(a)和隨密度變化(φ=5%)(b)的分布

3.3 剩余伽馬射線強度與地層密度的關系

設計脈沖中子源密度測井儀器時,其中一個最重要的參數就是源距L。伽馬探測器所放的位置,既要保證探測器的計數率較高又要保證反映地層密度的精度達到一定的需求。由3.2節(jié)中的模擬結果可知,在L≈55cm位置處,剩余伽馬射線強度不受孔隙度的影響,且對密度變化比較敏感,所以源距L最好選擇在50～60cm。

由3.1節(jié)中的模擬結果可知,在源距為50～60cm處,俘獲伽馬射線強度隨孔隙度的增大而減小,其變化關系如圖9a所示。由圖9a可見,俘獲伽馬射線強度與孔隙度φ呈指數關系。而由公式(3)可推知剩余伽馬射線強度J與孔隙度φ之間也呈指數關系。

選中源距為50～60cm之后,得到不同孔隙度條件下的剩余伽馬射線強度與密度的關系,其中L=55cm處的剩余伽馬射線強度與密度的關系曲線示于圖9b。使用模擬的數據進行擬合,得到剩余伽馬射線強度與地層密度、孔隙度以及源距的擬合關系為

(4)

其相關系數為0.97。需要說明的是,此關系式為每萬個中子產生的剩余伽馬射線強度,中子源能量為14MeV。在實際應用中,(4)式中的孔隙度φ為同一個儀器中子探測部分所獲得的孔隙度。

圖9 特定源距時俘獲伽馬射線強度Q與孔隙度φ的關系(a)和剩余伽馬射線強度J與密度ρ的關系(b)

3.4 剩余伽馬射線強度的影響因素

3.4.1 中子源強度N和能量E的影響

中子源的特性對剩余伽馬射線強度有著一定的影響,評價中子源的特性主要從源強度N和源能量E兩個方面考慮。圖10a是中子源強度對剩余伽馬射線強度影響(E=14MeV,φ=5%,VCl=10%,ρ=2.4g/cm3)的模擬結果,可以看出剩余伽馬射線強度隨中子源強度N的增大而增大。這是因為當中子源的強度增大時,與核素反應的中子數增加,產生的俘獲伽馬射線數量增大,從而剩余伽馬射線強度增加。圖10b是中子源能量對剩余伽馬射線強度影響(N=10000,φ=5%,VCl=10%,ρ=2.4g/cm3)的模擬結果,可以看出中子源的能量E對剩余伽馬射線強度幾乎沒有影響。

圖10 中子源強度(a)、中子源能量(b)和地層水含氯量(c)對剩余伽馬射線強度影響的模擬分析

3.4.2 地層水含氯量(VCl)的影響

氯是最重要的熱中子俘獲物質,俘獲截面大且放出的伽馬射線數量多,所以地層水中氯的含量對剩余伽馬射線強度有一定的影響,但影響不大。圖10c 為地層水含氯量對剩余伽馬射線影響(E=14MeV,N=10000,φ=5%,ρ=2.4g/cm3)的模擬結果。由圖10c可以看出,在源距選定之后,隨著水中含氯量的增大,剩余伽馬射線的強度也增大。這是因為含氯量增大使得中子與核素作用總的宏觀俘獲截面增加,產生的俘獲伽馬射線強度增大,從而剩余伽馬射線強度也隨之增大。

4 結論與認識

1) 脈沖中子源發(fā)射的快中子經地層減速作用產生的熱中子,以中子源為中心呈圓形分布；當地層孔隙度增大時,熱中子分布半徑變小。熱中子擴散作用導致熱中子云的半徑有所增大,而且隨著孔隙度的增加,熱中子的擴散長度減小。

2) 俘獲伽馬射線強度隨源距L的增大呈指數迅速降低。且當L>100cm時,俘獲伽馬射線強度幾乎為零；當L≈35cm時,俘獲伽馬射線強度與孔隙度無關；當L<35cm時,孔隙度大的地層俘獲伽馬射線強度大；當L>35cm時,孔隙度小的地層俘獲伽馬射線強度大。

3) 俘獲伽馬射線經過地層吸收之后,得到的剩余伽馬射線強度受孔隙度和地層密度的共同影響,隨密度的增加而減??；合適的源距為50～60cm,建立了剩余伽馬射線強度與地層密度、孔隙度以及源距的擬合關系,可以指導脈沖中子源密度測井儀器的設計。

4) 中子源強度、地層水含氯量對剩余伽馬射線強度也有影響,剩余伽馬射線強度隨中子源強度和地層水含氯量的增大而增大；而中子源能量對剩余伽馬射線強度幾乎沒有影響。

參 考 文 獻

[1] 于華偉,孫建孟,朱文娟,等.隨鉆脈沖中子密度測井的蒙特卡羅模擬研究[J].測井技術,2009,33(6):521-524

Yu H W,Sun J M,Zhu W J,et al.Monte-Carlo simulation method in pulsed neutron density logging while drilling[J].Well Logging Technology,2009,33(6):521-524

[2] Reichel N,Evans M,Allioli F,et al.Neutron-gamma density(ngd):principles,field test results and log quality control of a radioisotope-free bulk density measurement[J].SPWLA 53rdAnnual Logging Symposium, 2012,GGG

[3] Quirein J A,Jr H S,Chen D D,et al.Formation density prediction using pulsed neutron capture tools[J].SPWLA 46thAnnual Logging Symposium,2005,QQ

[4] Atfeh M,Daghar K A A,Marzouqi K A,et al.Neutron porosity and formation density acquisition without chemical sources in large carbonate reservoirs in the middle east-a case study[J].SPWLA 54thAnnual Logging Symposium,2013,KKK

[5] 吳文圣.Monte Carlo方法在核測井中的新應用[J].測井技術,2001,25(6):412-416

Wu W S.The uses of Monte Carlo method in nuclear logging[J].Well Logging Technology,2001,25(6):412-416

[6] 李舟波.鉆井地球物理勘探 [M].北京：地質出版社,2006：114-175

Li Z B.Geophysical exploration drilling [M].Beijing:Geological Publishing Press,2006：114-175

[7] 龐巨豐.核測井物理基礎[M].北京：石油工業(yè)出版社,2005：94-129

Pang J F.Physical basis of nuclear logging [M].Beijing:Petroleum Industry Press,2005:94-129

[8] 于華偉.隨鉆環(huán)境下脈沖中子測量地層密度的理論基礎研究[D].青島：中國石油大學(華東),2011

Yu H W.The fundamental research of the pulsed-neutron density logging while drilling[D].Qingdao:China University of Petroleum,2011

[9] Odom R C,Bailey S M,Wilson R D,et al.Pulsed neutron density measurements:modeling the depth of investigation and cased-hole wellborge uncertainties[J].SPWLA 40thAnnual Logging Symposium,1999,JJ