李國梁,張欽中,李雨蓮,張瓊,陳紅喜,劉曉斌

(1.中海油田服務股份有限公司油田技術(shù)研究院,天津300459;2.電子科技大學自動化工程學院,四川成都611731)

0 引 言

20世紀80年代以來,核測井數(shù)值模擬技術(shù)在國內(nèi)外逐步發(fā)展成熟,并在核測井原理研究、儀器設(shè)計、測量數(shù)據(jù)處理和標定、測井結(jié)果解釋等領(lǐng)域中發(fā)揮了重大作用[1-2]。蒙特卡羅方法被認為是解決復雜輻射輸運問題最精確的方法之一,具有較高的精度,被廣泛應用于石油勘探領(lǐng)域[3]。

核測井是根據(jù)井內(nèi)介質(zhì)、巖石及其孔隙流體的核物理性質(zhì),尋找油氣儲層的一種地球物理方法。中子孔隙度儀器將中子射入地層,根據(jù)近遠探測器計數(shù)率可以確定地層孔隙度值。

本文根據(jù)刻度井參數(shù)和儀器設(shè)計圖紙,使用蒙特卡羅軟件Geant4中的各種幾何體及布爾運算,完成了儀器和刻度井的物理模型構(gòu)建。并選用靈敏探測器以及合適的源和中子輸運物理過程,建立了隨鉆孔隙度測井的仿真模型。通過對比刻度井中儀器實測計數(shù)比值與仿真模型模擬計數(shù)比值,驗證了蒙特卡羅仿真模型的準確性。并運用平行幾何體加速蒙特卡羅模擬,降低了探測器計數(shù)比值的統(tǒng)計誤差,有效提高了模擬效率。

1 中子孔隙度測井原理

中子孔隙度測井是核測井中最常用的測量方法之一。中子通過與地層中的原子核發(fā)生彈性散射和非彈性散射損失能量,被撞粒子的質(zhì)量與中子越接近,則中子損失的能量越多。在地層元素中,氫與中子的能量最接近,是地層慢化中子的主要元素。因此,地層的含氫量決定了熱中子的空間分布[4]。而地層中的氫主要存在于地層孔隙中的油、氣、水中,因此,通過分析地層中氫對中子的影響程度,即可確定地層孔隙度。

根據(jù)熱中子擴散理論,得到近遠源距探測器熱中子計數(shù)比為[5]

(1)

式中,Nn、Nf分別為近、遠源距探測器熱中子計數(shù)率,cps;Le、Lt分別為快中子與熱中子的減速長度和擴散長度,mm;r1、r2為近、遠探測器的源距,mm;R為近、遠源距探測器熱中子計數(shù)比值。由式(1)可見,探測器計數(shù)比值只與地層減速性質(zhì)有關(guān)。

熱中子通過測井儀器上的電子線路系統(tǒng)產(chǎn)生的脈沖進行計數(shù),因此,從記錄的脈沖中可以得到地層孔隙度的相關(guān)信息[6]。

2 蒙特卡羅模擬及驗證

2.1 蒙特卡羅模型的建立

Geant4是一款由歐洲核子中心主導,世界范圍內(nèi)的相關(guān)核與粒子物理學家參與開發(fā)和維護的模擬粒子與物質(zhì)相互作用的蒙特卡羅應用計算程序的開源工具包[7-8]。其模擬程序主要包括:幾何模型構(gòu)建、幾何材料定義、初級以及次級粒子描述、初始事件的生成、粒子跟蹤、粒子與材料相互作用的物理過程描述、事件數(shù)據(jù)的生成、事件數(shù)據(jù)的保存、粒子徑跡的可視化等。目前已經(jīng)廣泛應用于高能物理、原子核實驗物理、探測器模擬、核醫(yī)學、航空航天、加速器物理、輻射防護等領(lǐng)域[8]。

INP675和INP800儀器是中海油田服務股份有限公司自主研發(fā)的2種隨鉆中子測井儀器,INP675儀器用于8.5 in(1)非法定計量單位,1 in=2.54 cm,下同的井眼,INP800用于12.25 in的井眼。測井時INP675和INP800中子測井儀器均置于鉆鋌中,可以實現(xiàn)地層孔隙度的隨鉆測量。放射源使用Am-Be同位素中子源。INP675儀器的2組探測器內(nèi)均有6根成圓弧形排列的3He管,INP800儀器探測器的結(jié)構(gòu)與INP675相似,每組探測器內(nèi)有3根成圓弧形排列3He管。測量環(huán)境為中海油田服務股份有限公司的刻度井群,包括巖石、鋼桶、鼠洞等結(jié)構(gòu)。根據(jù)中海油田服務股份有限公司《中子孔隙度刻度井群測試規(guī)范》[9]完成實測與建模。

仿真模型主要包括井基座、鋼桶、井體、鉆鋌等部分,源和3He管支架處于鉆鋌內(nèi)部。鉆鋌中心留有水孔,用于輸送鉆井液。基座大致分為3部分:井體、鼠洞、實心底座。在Geant4中提供了很多現(xiàn)成的幾何形狀,如球體、長方體、錐體等。基座建模時,使用軟件內(nèi)部的庫函數(shù)構(gòu)建相應大小的圓柱體,并運用布爾減運算做出鼠洞。確定幾何體形狀后,設(shè)置基座的材料等屬性,完成邏輯體的建立。最后,確定基座的母體和擺放位置,將邏輯體轉(zhuǎn)化為物理體,完成基座建模。其余物理體模型也按照上述步驟,參照參考文獻[9]和儀器設(shè)計圖紙,使用Geant4庫函數(shù)中提供的幾何體形狀和3種布爾運算(交、并、補)完成建模。

儀器的源為Geant4工具包中的GPS(General Particle Source)類定義的各向同性的標準Am-Be化學源。選擇屏蔽物理過程模擬井下物質(zhì)對中子的慢化和吸收作用。使用靈敏探測器記錄接收到的中子計數(shù)。將探測器探測的粒子沉積能量記錄在energyDeposit變量中,當energyDeposit大于0時,表明3He管探測器中有熱中子進入,探測器計數(shù)值加1,以此實現(xiàn)探測器的計數(shù)。

INP675與INP800儀器結(jié)構(gòu)和建模方法相似,僅在工具大小、源距和3He管數(shù)量等儀器參數(shù)有些許差距。以R7刻度井和INP800儀器為例,井體和儀器的結(jié)構(gòu)如圖1所示。

圖1 刻度井現(xiàn)場及仿真模型示意圖

2.2 蒙特卡羅模型驗證

選取了20口不同的灰?guī)r刻度井進行模擬,其孔隙度在9.3～50.8 p.u.區(qū)間內(nèi)變化,得到每種儀器在刻度井中的模擬數(shù)據(jù),記錄近遠源距探測器計數(shù),并計算近遠源距探測器熱中子計數(shù)比值R,通過與儀器實測計數(shù)比值的對比,驗證仿真模型的可靠性和準確性。部分刻度井信息見表1,表1中孔隙度不確定度代表刻度井的實際孔隙度在孔隙度標稱值附近的變化量,如S5井的實際孔隙度為9.3±0.1 p.u.。

表1 部分灰?guī)r刻度井信息

分別在10口刻度井中,直接對2種儀器的實測計數(shù)比值和模擬計數(shù)比值進行匹配,得到的線性相關(guān)性見圖2。

圖2 儀器實測計數(shù)比值和模擬計數(shù)比值擬合結(jié)果對比圖

根據(jù)擬合結(jié)果,刻度井中模擬計數(shù)比值和實測計數(shù)比值的相關(guān)性分別為0.974 6和0.952 5,模擬值和實測值有較強的線性擬合關(guān)系,證明仿真模型能較好地反映實際的井下測量環(huán)境。為了進一步評價刻度井中模擬數(shù)據(jù)和實測數(shù)據(jù)的差異,使用最小二乘法,對數(shù)據(jù)進行線性擬合,擬合結(jié)果見圖3。

圖3 儀器實測計數(shù)比值和模擬計數(shù)比值對比

儀器模擬計數(shù)比值和實測計數(shù)比值具有較好的線性關(guān)系。經(jīng)過多次實驗,當模擬源出射粒子數(shù)為5×108時,探測器計數(shù)比值的統(tǒng)計誤差小于1%,可以達到較高的精度要求。計算儀器實測計數(shù)比值和模擬計數(shù)比值的相對誤差,見表2和圖4。

圖4 儀器實測計數(shù)比值和模擬計數(shù)比值的相對誤差

表2 儀器實測計數(shù)比值和模擬計數(shù)比值的相對誤差

不同刻度井中INP675儀器的實測計數(shù)比值和模擬計數(shù)比值的最大相對誤差為4.48%,最小為2.31%,平均相對誤差為3.46%;INP800儀器計數(shù)比值的實測計數(shù)比值和模擬計數(shù)比值的相對誤差最大為2.15%,最小為1.61%,平均相對誤差為1.91%,模擬值與實測值較為接近。綜合2種儀器的刻度井實測計數(shù)比值和模擬計數(shù)比值的對比結(jié)果,2種儀器實測值和模擬值的相關(guān)性均大于0.950 0,相對誤差均低于4.50%,到達了較好的匹配效果。

以上數(shù)據(jù)分析充分說明仿真模型的模擬值和實測值有著較好的匹配度,證明仿真模型中井下環(huán)境的建模、儀器的構(gòu)建、物理過程的選擇等是恰當?shù)?。通過對不同環(huán)境參數(shù)條件下的仿真模擬,可以推算出復雜井況下探測器響應與孔隙度的關(guān)系,從而節(jié)省大量的人力物力。

3 蒙特卡羅方法優(yōu)化研究

由于探測器體積相對于整個地層較小,同時源是各向同性的,大部分粒子不能被探測器探測到,因此,傳統(tǒng)的蒙特卡羅方法需要大量的粒子和時間進行模擬才能達到預期的精度水平。以INP800儀器為例展開研究,該儀器在灰?guī)r孔隙度為50.8%時,且當源出射粒子數(shù)為108時,計數(shù)比值的統(tǒng)計誤差小于5%,屬于可接受的誤差范圍。經(jīng)實驗,使用Inter?xeon(R)Platinum 8280L處理器的服務器單線程運行需要32 d。因此,加速蒙特卡羅模擬對于核測井至關(guān)重要。

在粒子輸運過程中,減方差方法是提高仿真效率的方法之一。本文在仿真中使用的減方差方法是基于Geant4的平行幾何體實現(xiàn)的[10],通過創(chuàng)建與物理幾何體相關(guān)聯(lián)的平行幾何體,引入新的概率密度函數(shù)對粒子使用分裂或輪盤賭技巧,實現(xiàn)偏倚仿真,達到減方差的目的[11]。

以本設(shè)計中使用的平行幾何體為例(見圖5)。對幾何體進行20次切割,形成了多個具有不同重要性的區(qū)域。重要性I代表的是對該區(qū)域的感興趣程度,也就是該區(qū)域粒子對探測器計數(shù)的貢獻情況。

圖5 平行幾何體示意圖

(2)

式中,Qc為該區(qū)域內(nèi)被探測器探測到粒子的總權(quán)重;Wt為該區(qū)域所有粒子的總權(quán)重。因為在中子孔隙度測井模型中希望增強從地層中散射回儀器的中子數(shù)量,所以儀器附近區(qū)域的重要性會高于遠離儀器的區(qū)域的重要性。根據(jù)各區(qū)域重要性大小定義生存概率P

(3)

式中,Ipost、Ipre分別為粒子軌跡終點方向和起點方向區(qū)域的重要性。粒子經(jīng)過重要性不同的2個區(qū)域時:當P<1時,粒子分裂;當P>1時,由輪盤賭決定粒子被殺死還是存活。無論發(fā)生粒子分裂還是粒子存活,其權(quán)重W都會發(fā)生改變,變?yōu)閃·P。

原仿真與偏倚仿真的對比結(jié)果見表3,其中運行時間為單線程運行時間,計數(shù)比值的相對統(tǒng)計誤差(σRE)計算方法見式(4)[11]。在偏倚仿真中,由于粒子權(quán)重發(fā)生了變化,計算相對統(tǒng)計誤差時應該使用有效樣本數(shù)[12][見式(5)],而不能使用無權(quán)重樣本數(shù)。在蒙特卡羅仿真中一般使用品質(zhì)因子FOM[見式(6)]衡量仿真效率,它同時考慮了相對統(tǒng)計誤差與運算時間,數(shù)值越大代表效率越高。

表3 粒子數(shù)為108時裸眼仿真輸出

(4)

(5)

(6)

式中,σRE為近、遠源距探測器計數(shù)比值的相對統(tǒng)計誤差,%;en、ef分別為近、遠源距探測器計數(shù)的相對統(tǒng)計誤差,%;Nn、Nf分別為近、遠源距探測器計數(shù);Wi為第i個粒子的權(quán)重;Neff為有效樣本數(shù);t為單個仿真模型的運算時間,min。

由表3可見,應用平行幾何體后,儀器計數(shù)比值的相對統(tǒng)計誤差是原來的1/5;雖然偏倚仿真的單線程運行時間比原仿真長,但它的相對統(tǒng)計誤差大大降低。因此,仿真效率的評價應使用品質(zhì)因子,綜合考慮相對統(tǒng)計誤差和運行時間這2種因素。偏倚仿真的品質(zhì)因子至少高于原仿真結(jié)果的50倍,說明采用偏倚技術(shù)可有效提高仿真效率。圖6為地層中的粒子分布圖且數(shù)據(jù)都已標準化,圖6中顏色越接近黃色,表示粒子密度越高,反之,顏色越接近藍色,表示粒子密度越低。結(jié)果表明,采用平行幾何體時,粒子分布更廣,可探測到的粒子更多。這一方法的應用可以提高后續(xù)大規(guī)模模型運算的效率,節(jié)約計算成本。

圖6 偏倚仿真模擬與原仿真模擬x-y軸粒子分布俯視圖

4 結(jié) 論

(1)使用Geant4構(gòu)建了INP675和INP800這2種隨鉆中子測井儀和20口刻度井的仿真模型,并選用靈敏探測器和合適的物理過程建立了隨鉆中子儀器孔隙度測量的仿真模型。

(2)刻度井中,INP675儀器和INP800儀器的模擬計數(shù)比值和實測計數(shù)比值的相關(guān)性分別為0.974 6和0.952 5,相對誤差均小于4.5%,達到了較好的匹配結(jié)果,驗證了仿真模型的準確性和可靠性。

(3)應用平行幾何體加速蒙特卡羅模擬,儀器計數(shù)比值的相對統(tǒng)計誤差是原模擬結(jié)果的1/5,品質(zhì)因子至少高于原仿真結(jié)果的50倍,有效提高了仿真效率,節(jié)約了成本。

(4)套管井會對中子輸運過程和中子儀器測量產(chǎn)生更為復雜的影響。完成了INP675、INP800儀器在裸眼井中孔隙度測量的準確模擬,為過套管孔隙度測量的環(huán)境校正及相關(guān)反演算法的研究提供了基礎(chǔ)。