基于卷積神經(jīng)網(wǎng)絡(luò)和MPGA-LM算法的陣列側(cè)向測(cè)井快速反演方法

吳易智， 范宜仁*， 巫振觀， 鄧少貴，張盼， 陳詩(shī)宇， 尹中旭

1 中國(guó)石油大學(xué)(華東)地球科學(xué)與技術(shù)學(xué)院， 青島 266580 2 海洋國(guó)家實(shí)驗(yàn)室海洋礦產(chǎn)資源評(píng)價(jià)與探測(cè)技術(shù)功能實(shí)驗(yàn)室， 青島 266237

0 引言

憑借著較高的分辨率以及豐富的探測(cè)信息，陣列側(cè)向測(cè)井被廣泛用于大斜度井/水平井的油氣藏評(píng)價(jià)(Rasmus, 1982; Itskovich et al., 1998; 鄧少貴等, 2010).然而在致密性儲(chǔ)層電阻率評(píng)價(jià)過(guò)程中，由于地層裂縫和薄互層發(fā)育以及沉積顆粒的大小、磨圓不一等因素，導(dǎo)致儲(chǔ)層各向異性特征明顯，測(cè)井響應(yīng)復(fù)雜，為油氣層的識(shí)別與評(píng)價(jià)帶來(lái)一定的難度(Salazar and Torres-Verdín, 2009; 高杰等, 2010; Gao et al., 2013).因此，如何從測(cè)井響應(yīng)中準(zhǔn)確提取地層真實(shí)電阻率成為電測(cè)井領(lǐng)域亟待解決的難題之一.

非線性反演技術(shù)可以充分考慮井眼、侵入等因素對(duì)測(cè)井響應(yīng)的影響，成為大斜度井/水平井井旁電阻率剖面重構(gòu)的有效手段(Nam et al., 2010; Wang, 2011; 潘克家等, 2013).其算法包括確定性算法和隨機(jī)類算法.其中，確定性算法主要基于最小二乘理論，其核心在于通過(guò)不斷迭代模型參數(shù)的方式使得模擬數(shù)據(jù)與觀測(cè)數(shù)據(jù)的偏差最小(Osman and Laporte, 1996; Habashy and Abubakar, 2004; Madsen et al., 2004; 姚東華等, 2010; Hu and Fan, 2018)，該類算法包括最速下降法、高斯牛頓法以及列文伯格馬奎特法等(Anderson, 2001).由于能夠有效消除圍巖、層厚和侵入等因素的影響，確定性算法被廣泛用于2D陣列側(cè)向測(cè)井資料的快速處理(Wang et al., 2009; 王磊等, 2018).然而對(duì)于斜井各向異性等三維地層環(huán)境而言，正演效率過(guò)低，導(dǎo)致反演迭代過(guò)程中雅可比矩陣或海森矩陣求取耗時(shí)，同時(shí)該類算法對(duì)初值依賴性強(qiáng)，易陷于局部極小值，不利于陣列側(cè)向測(cè)井資料的快速精確處理(Pardo et al., 2007; Dubourg et al., 2017; Wu et al., 2020).不同于確定性算法，隨機(jī)類算法由單點(diǎn)或多點(diǎn)啟發(fā)，搜索步長(zhǎng)和方向具有隨機(jī)性，可以在求解域內(nèi)全空間搜索，全局尋優(yōu)能力強(qiáng)，進(jìn)而能夠提高資料處理精度，此類算法包括粒子群優(yōu)化算法、模擬退火算法以及差分進(jìn)化算法等(Eiben et al., 1994;Yang and Torres-Verdín, 2011; Li et al., 2020).然而，兩類算法在尋優(yōu)過(guò)程中皆要不斷調(diào)用正演，尤其對(duì)于三維地層環(huán)境而言，正演速度仍是制約反演效率的關(guān)鍵因素.

近年來(lái)，隨著智能算法的興起，有研究學(xué)者將深度學(xué)習(xí)技術(shù)引入地球物理反演當(dāng)中(林年添等, 2018; Shahriari et al., 2020; Wu and Fan, 2021; Zhang et al.,2021)，以期能夠解決反演策略中正演效率低的難題.其中，反向傳播神經(jīng)網(wǎng)絡(luò)(BPNN)作為一種強(qiáng)大的監(jiān)督類學(xué)習(xí)架構(gòu)，自從提出以來(lái)被廣泛用于特征識(shí)別以及函數(shù)逼近領(lǐng)域(Rumelhart and McClelland, 1986; Xu et al., 2018).然而，面對(duì)復(fù)雜輸入輸出問(wèn)題，BPNN存在超參數(shù)優(yōu)化困難、遷移性差、訓(xùn)練參數(shù)過(guò)多以及不易收斂等缺點(diǎn).為此，深層卷積神經(jīng)網(wǎng)絡(luò)(CNN)憑借著局部連接和權(quán)重共享的優(yōu)點(diǎn)逐漸代替BPNN(Fukushima, 1980; LeCun et al., 1998; Krizhevsky et al., 2017).在測(cè)井領(lǐng)域，Zhu等(2018)通過(guò)建立五層卷積神經(jīng)網(wǎng)絡(luò)實(shí)現(xiàn)了地層巖性特征的快速識(shí)別；Li等(2019)基于1D模型和多層卷積神經(jīng)網(wǎng)絡(luò)實(shí)現(xiàn)了大斜度井/水平井隨鉆電磁波測(cè)井的快速反演；Jin等(2019)通過(guò)優(yōu)化卷積神經(jīng)網(wǎng)絡(luò)的損失函數(shù)實(shí)現(xiàn)了隨鉆測(cè)井領(lǐng)域的實(shí)時(shí)地質(zhì)導(dǎo)向.然而，上述研究皆是針對(duì)簡(jiǎn)單輸入輸出問(wèn)題，而對(duì)于陣列側(cè)向測(cè)井3D反演問(wèn)題，國(guó)內(nèi)外尚無(wú)研究.因此，針對(duì)陣列側(cè)向測(cè)井資料處理，如何提升高維正演效率，如何保證反演精度是當(dāng)前智能測(cè)井領(lǐng)域亟需解決的難題.

為此，本文從三維有限元技術(shù)(3D-FEM)出發(fā)，首先對(duì)斜井各向異性地層參數(shù)進(jìn)行敏感性分析，確定侵入、各向異性和地層傾角對(duì)陣列側(cè)向測(cè)井的敏感性大??；其次，建立三維陣列側(cè)向測(cè)井響應(yīng)數(shù)據(jù)庫(kù)，引入二維卷積神經(jīng)網(wǎng)絡(luò)(2D-CNN)技術(shù)，并基于二進(jìn)制轉(zhuǎn)換將地層模型轉(zhuǎn)為網(wǎng)絡(luò)輸入層所需的2D圖像，提高網(wǎng)絡(luò)識(shí)別局部特征的能力；接著，建立深層2D-CNN結(jié)構(gòu)，并引入“丟棄”(Dropout)和“填充”(Padding)機(jī)制，實(shí)現(xiàn)三維陣列側(cè)向快速正演；最后，基于敏感性分析和快速正演，將多種群遺傳算法與列文伯格馬奎特算法相結(jié)合(MPGA-LM)，對(duì)傾斜各向異性地層的電阻率剖面進(jìn)行了快速精確重構(gòu).

1 陣列側(cè)向測(cè)井?dāng)?shù)值模擬

1.1 三維有限元數(shù)值模擬

本文以斯倫貝謝研發(fā)的高分辨率陣列側(cè)向測(cè)井儀HRLA為例，儀器包括一個(gè)主電極A0、六對(duì)屏蔽電極A1-A6(A1′-A6′)及兩對(duì)監(jiān)督電極M1、M2(M1′、M2′)，屏蔽電極和監(jiān)督電極關(guān)于主電極對(duì)稱，通過(guò)改變不同電極間的收發(fā)組合方式，可以形成6種具有不同探測(cè)深度的視電阻率曲線，探測(cè)深度由淺及深分別命名為RLA0，RLA1，RLA2，RLA3，RLA4，RLA5.

不同探測(cè)模式視電阻率可以通過(guò)測(cè)量監(jiān)督電極電位進(jìn)行計(jì)算，其中淺探測(cè)模式RLA0需要測(cè)量同側(cè)兩個(gè)監(jiān)督電極的電勢(shì)差值ΔUM1M2(RLA0)，其他探測(cè)模式只需測(cè)量任意監(jiān)督電極上電勢(shì)，其視電阻率計(jì)算公式如下：

(1)

(2)

式中，i為1～5，分別代表上述5種工作方式；Ki為第i種測(cè)井模式的電極系數(shù)，可在均勻地層條件下獲得；RRLAi為第i種測(cè)井模式的視電阻率；UM1(RLAi)為第i種測(cè)井模式的主電極電位；I0(RLAi)為第i種測(cè)井模式的主電極電流.

可以看出，求解陣列側(cè)向測(cè)井視電阻率的關(guān)鍵在于得到監(jiān)督電極上的電勢(shì)，為此必須將整個(gè)空間的電位分布函數(shù)求出.由于陣列側(cè)向測(cè)井用的是低頻交流電，故可視為穩(wěn)定電流，進(jìn)而該問(wèn)題轉(zhuǎn)為穩(wěn)流場(chǎng)電勢(shì)分布函數(shù)的求解問(wèn)題.

根據(jù)邊界條件和變分原理，陣列側(cè)向測(cè)井響應(yīng)可歸結(jié)為求取泛函極值問(wèn)題，泛函構(gòu)造式如下：

(3)

式中，Λ為除去電極部分的求解區(qū)域，Ii和Ui分別代表儀器每個(gè)電極的電流和電勢(shì)大小.

在電阻率不同的子空間交界面上，電流、電位具有連續(xù)性，滿足如下邊界條件：

(4)

式中，“-”表示交界面左側(cè)(上側(cè))，“+”表示交界面右側(cè)(下側(cè))；n為交界面單位法向量.

在無(wú)窮遠(yuǎn)邊界和絕緣電極表面還需滿足Dirichlet和Neumann邊界條件，如下所示：

U|Γ1=0,

(5)

式中，Γ1和Γ2分別為無(wú)窮遠(yuǎn)邊界和絕緣電極表面.

在主電極和屏蔽電極表面，滿足等值面邊界條件：

(6)

式中，UAi為電極表面電勢(shì)，Γ3主電極或屏蔽電極表面，n為電極表面單位法向量，Ii為電極發(fā)射電流，i為0～5.

為求解上述泛函的最小值，一般采用三維有限元方法(3D-FEM)，最后采用前線解法求取含有大型稀疏矩陣的線性方程組，如式(7)，該方法由于可以在求解過(guò)程中對(duì)剛度矩陣邊消元邊安裝，可有效提升正演計(jì)算速度(張庚驥, 2009).

KΦ=F,

(7)

1.2 影響因素分析

陣列側(cè)向測(cè)井受泥漿侵入、地層傾角以及各向異性等多種因素影響，因此，確定各因素對(duì)測(cè)井響應(yīng)影響強(qiáng)弱對(duì)于反演初值選取、約束條件施加具有重要意義.本節(jié)系統(tǒng)分析了侵入深度、各向異性和地層相對(duì)傾角等因素對(duì)陣列側(cè)向測(cè)井響應(yīng)的敏感性，定義敏感因子S如下：

(8)

其中，i為陣列側(cè)向探測(cè)模式，此處取1～5，x為地層參數(shù)，分別代表侵入深度Di，各向異性系數(shù)λ和地層傾角θ.

首先，分析侵入對(duì)陣列側(cè)向測(cè)井響應(yīng)的影響.以無(wú)限厚一層模型為例，地層為各向同性，侵入帶電阻率為2 Ωm，侵入深度Di變化范圍為0～3 m，侵入帶與原狀地層電阻率對(duì)比度范圍Rat為1～20.敏感性分布圖如圖1所示，橫縱坐標(biāo)分別為侵入深度和電阻率對(duì)比度.可以看出，陣列側(cè)向測(cè)井對(duì)侵入深度敏感因子最大約為2，且隨著探測(cè)深度的增大，對(duì)侵入深度的敏感區(qū)間也隨之增大：淺、中、深探測(cè)模式RLA1、RLA3和RLA5的敏感區(qū)域分別為0.2～0.5 m、0.3～1.0 m和0.6～2.3 m.

圖1 HRLA對(duì)侵入深度的敏感性 (a)?RLA1/?Di； (b)?RLA3/?Di； (c)?RLA5/?DiFig.1 The sensitivity of HRLA to invasion depth

接著，分析各向異性和地層傾角對(duì)陣列側(cè)向測(cè)井的影響.以無(wú)限厚模型為例，地層為各向異性，地層水平電阻率為2 Ωm，各向異性系數(shù)變化范圍為1～4，地層相對(duì)傾角變化范圍為0°～90°.圖2和圖3分別表示各向異性和地層相對(duì)傾角的敏感性分布圖，橫縱坐標(biāo)分別為各向異性系數(shù)和地層相對(duì)傾角.可以看出：(1)陣列側(cè)向測(cè)井對(duì)各向異性系數(shù)和地層相對(duì)傾角的敏感因子最大約為1；(2)陣列側(cè)向測(cè)井在高角度地層(60°～90°)對(duì)各向異性最為敏感，且探測(cè)深度較深的測(cè)井模式對(duì)各向異性的敏感性較大，原因在于當(dāng)?shù)貙酉鄬?duì)傾角逐漸增大時(shí)，儀器電極發(fā)射電流方向逐漸朝著垂直電阻率方向傾斜.因此，垂直電阻率對(duì)測(cè)井響應(yīng)貢獻(xiàn)變大，各向異性敏感性增強(qiáng)，而深探測(cè)曲線因探測(cè)范圍最廣，受其影響最大；(3)陣列側(cè)向測(cè)井在中高角度(50°～80°)對(duì)地層相對(duì)傾角最為敏感，且不同測(cè)井模式對(duì)相對(duì)傾角的敏感性與各向異性相似，皆是探測(cè)深度較深的測(cè)井模式對(duì)其敏感性較大，原因在于地層傾角對(duì)陣列側(cè)向曲線的影響本質(zhì)是因?yàn)榈貙痈飨虍愋缘挠绊?，因?yàn)榈貙觾A角增大，導(dǎo)致垂直電阻率的影響增大，進(jìn)而使得陣列側(cè)向五條曲線出現(xiàn)變化.

圖2 HRLA對(duì)各向異性的敏感性 (a)?RLA1/? λ； (b)?RLA3/? λ； (c)?RLA5/? λ.Fig.2 The sensitivity of HRLA to anisotropy

圖3 HRLA對(duì)地層相對(duì)傾角的敏感性 (a)?RLA1/? θ； (b)?RLA3/? θ； (c) ?RLA5/? θ.Fig.3 The sensitivity of HRLA to dipping angle

2 深層卷積神經(jīng)網(wǎng)絡(luò)快速正演模型

2.1 深層卷積神經(jīng)網(wǎng)絡(luò)理論

深層卷積神經(jīng)網(wǎng)絡(luò)憑借著強(qiáng)大的泛化特性以及學(xué)習(xí)能力被廣泛用于圖像分類、語(yǔ)音檢測(cè)以及物體識(shí)別等復(fù)雜特征識(shí)別領(lǐng)域.一般而言，卷積神經(jīng)網(wǎng)絡(luò)包括輸入層、卷積層、匯聚層(池化層)、全連接層和輸出層，本文以經(jīng)典二維卷積神經(jīng)網(wǎng)絡(luò)——LeNet-5為例，其結(jié)構(gòu)如圖4所示.其中，輸入層為三通道RGB圖像或者單通道灰度圖像，卷積層包含多個(gè)卷積核，不同的卷積核相當(dāng)于不同的特征提取器，結(jié)構(gòu)為一個(gè)二維矩陣，圖像經(jīng)過(guò)卷積層可提取多個(gè)局部特征并保存；匯聚層的作用在于特征選擇，降低特征數(shù)量，從而減少參數(shù)數(shù)量，其核心思想在于步長(zhǎng)間隔的選擇，達(dá)到“降維”的作用；全連接層與全連接前饋神經(jīng)網(wǎng)絡(luò)相同，不再贅述.

圖4 LeNet-5網(wǎng)絡(luò)結(jié)構(gòu)Fig.4 The architecture of LeNet-5

相較于全連接神經(jīng)網(wǎng)絡(luò)，卷積神經(jīng)網(wǎng)絡(luò)最大的特點(diǎn)在于局部連接和權(quán)重共享，如圖5所示，圖5a為傳統(tǒng)的二層全連接神經(jīng)網(wǎng)絡(luò)，包括3個(gè)輸入，5個(gè)輸出，需訓(xùn)練參數(shù)包括15個(gè)權(quán)重和1個(gè)偏置共16個(gè)參數(shù)；圖5b為二層卷積神經(jīng)網(wǎng)絡(luò)，包括3個(gè)輸入，5個(gè)輸出，由于輸入層每個(gè)神經(jīng)元只需要與下層局部窗口內(nèi)的神經(jīng)元連接，且每個(gè)窗口內(nèi)權(quán)重共享，若卷積核大小為3，則需訓(xùn)練參數(shù)包括3個(gè)權(quán)重和1個(gè)偏置共4個(gè)參數(shù)，可大大減少訓(xùn)練參數(shù)的數(shù)量，提高訓(xùn)練效率.

圖5 神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu) (a) 全連接層； (b) 卷積層.Fig.5 The connection feature of network (a) Fully connected layer； (b) Convolution layer.

同時(shí)，為提高深層神經(jīng)網(wǎng)絡(luò)的訓(xùn)練效率，通常卷積層和全連接層分別引入填充(Padding)和丟棄(Dropout)機(jī)制，如圖6所示，通過(guò)在圖像卷積之前填充一層像素，使卷積之后的圖片與原圖尺寸相同，進(jìn)而使得卷積層得以加深，提高訓(xùn)練精度，如圖6a；另外通過(guò)隨機(jī)丟棄部分神經(jīng)元，可有效防止網(wǎng)絡(luò)過(guò)擬合，加快訓(xùn)練進(jìn)程，如圖6b.

圖6 神經(jīng)元連接機(jī)制 (a) “Padding”機(jī)制； (b) “Dropout”機(jī)制.Fig.6 Neuronal connection mechanism (a) “Padding” mechanism； (b) “Dropout” mechanism.

2.2 數(shù)據(jù)庫(kù)構(gòu)建及模型可視化

神經(jīng)網(wǎng)絡(luò)模型的可靠性和穩(wěn)定性很大程度上依賴于訓(xùn)練數(shù)據(jù)的精度，因此，本文基于嚴(yán)格三維有限元算法，構(gòu)建大型數(shù)據(jù)集，并將其載入深層2D-CNN進(jìn)行訓(xùn)練，以此建立快速回歸模型.考慮圍巖和層厚對(duì)陣列側(cè)向測(cè)井響應(yīng)的影響，數(shù)據(jù)集基于三層地層模型構(gòu)建，如圖7a所示，該模型亦作為2D-CNN的輸入.具體地，令兩側(cè)圍巖為各向同性地層，電阻率相同，中間層為有侵各向異性地層，模型參數(shù)范圍如下：各向異性系數(shù)λ為1～4，分為13組，間隔0.25；侵入深度Di為0～1 m，分為11組，間隔0.1 m；圍巖電阻率Rs、沖洗帶水平電阻率Rxoh和目的層水平電阻率Rth范圍皆為0.1～500 Ωm，將其以對(duì)數(shù)形式各分為20組，對(duì)數(shù)間隔為0.19；兩側(cè)圍巖為半無(wú)限厚地層；目的層厚度H為0～5 m，分為16組，間隔0.33 m；測(cè)量點(diǎn)為以目的層中點(diǎn)為中心，上下各取兩個(gè)點(diǎn)，共5個(gè)測(cè)量點(diǎn)，采樣點(diǎn)間隔為0.1 m.需要注意的是，為簡(jiǎn)化數(shù)據(jù)庫(kù)構(gòu)建流程，將其按地層傾角分為10個(gè)子庫(kù)，每個(gè)子庫(kù)地層傾角間隔10°.因此，每個(gè)子庫(kù)有13×11×203×16，共18,304,000個(gè)樣本，每個(gè)樣本包括六個(gè)模型參數(shù)，共109,824,000個(gè)模型數(shù)據(jù)，將其儲(chǔ)存至變量M中，大小為18304000×6；同時(shí)每個(gè)樣本含有五個(gè)測(cè)量點(diǎn)，每個(gè)測(cè)量點(diǎn)包括五個(gè)電阻率值(RLA1～RLA5)，因此，每組共包含457,600,000個(gè)數(shù)據(jù)點(diǎn)，將數(shù)據(jù)點(diǎn)作為輸出儲(chǔ)存至變量D中，大小為18304000×25.在計(jì)算機(jī)內(nèi)存中，每個(gè)數(shù)據(jù)點(diǎn)占8個(gè)字節(jié)，因此，每組數(shù)據(jù)大小約為0.944 GB，數(shù)據(jù)庫(kù)總大小為9.44 GB.

對(duì)于復(fù)雜模型回歸問(wèn)題，1D-CNN訓(xùn)練精度無(wú)法滿足需求，因此，本文采用學(xué)習(xí)和泛化能力較強(qiáng)的2D-CNN結(jié)構(gòu).而對(duì)于2D-CNN而言，輸入層普遍為2D圖像，為滿足網(wǎng)絡(luò)需求，本文基于二進(jìn)制轉(zhuǎn)換將一系列一維模型數(shù)組轉(zhuǎn)化為二維字符串，進(jìn)而合成一系列二值圖像(Zhong et al., 2019).具體轉(zhuǎn)換過(guò)程如圖7所示，其中圖7a為三層地層模型，具體模型參數(shù)如下：兩側(cè)圍巖為半無(wú)限厚各向同性地層，電阻率為2.564 Ωm，中間層層厚為3.156 m，侵入深度為0.428 m，侵入帶和目的層水平電阻率分別為11.278 Ωm和47.433 Ωm，各向異性系數(shù)為1.3，圖7b為模型參數(shù)對(duì)應(yīng)的二進(jìn)制矩陣，圖7c為基于二進(jìn)制矩陣轉(zhuǎn)換為的二值圖像，每個(gè)模型可唯一對(duì)應(yīng)一個(gè)二值圖像.由于數(shù)據(jù)庫(kù)中最大數(shù)值為8000，且小數(shù)后保留3位，故確定二進(jìn)制字符串長(zhǎng)度為32，若字符串長(zhǎng)度過(guò)長(zhǎng)，可自動(dòng)補(bǔ)零.

圖7 模型可視化過(guò)程 (a) 地層模型； (b) 模型參數(shù)對(duì)應(yīng)二進(jìn)制矩陣； (c) 二值圖像.Fig.7 The process of model visualization (a) Formation model； (b) The binary matrix correspond to model parameters； (c) Binary image corresponding to formation model.

2.3 2D-CNN損失函數(shù)及數(shù)據(jù)歸一化

由于前饋神經(jīng)網(wǎng)絡(luò)的局限性導(dǎo)致訓(xùn)練精度偏低，為此本文將反向傳播技術(shù)應(yīng)用于2D-CNN網(wǎng)絡(luò)架構(gòu)中，定義損失函數(shù)如下：

(9)

g(w,b,x)=w?x+b,

(10)

σ(z)=max(z,0),

(11)

(12)

式中，L為神經(jīng)網(wǎng)絡(luò)層數(shù)，Kl和Ml分別為第l層卷積核數(shù)量和卷積核尺寸.

接著，通過(guò)計(jì)算損失函數(shù)關(guān)于每層參數(shù)的梯度進(jìn)行反向傳播，為不失一般性，本文使用鏈?zhǔn)角髮?dǎo)法則對(duì)卷積層第l層的參數(shù)進(jìn)行求導(dǎo):

(13)

由式(10)和式(11)，易得

(14)

(15)

式中，k為第l層第k個(gè)特征映射，U為第l+1層的卷積核數(shù)量，X(l,k)=σ(g(l,k)).

同時(shí)，由于每組樣本中有5個(gè)測(cè)量點(diǎn)，即Y= (Y1，Y2，Y3，Y4，Y5)T，每個(gè)測(cè)量點(diǎn)包含五個(gè)視電阻率值，即Yi=(RLAi1,RLAi2,RLAi3,RLAi4,RLAi5)T.為提升訓(xùn)練精度，可將數(shù)據(jù)進(jìn)行歸一化處理，本文采用最小最大值歸一化方法，將輸出數(shù)據(jù)統(tǒng)一歸一化至[0,1]，公式如下：

(16)

需要注意的是，對(duì)于卷積神經(jīng)網(wǎng)絡(luò)模型預(yù)測(cè)的數(shù)據(jù)而言，需將其反歸一化至正常值.

2.4 深層卷積神經(jīng)網(wǎng)絡(luò)設(shè)計(jì)及測(cè)試

本文以9層卷積神經(jīng)網(wǎng)絡(luò)為例，對(duì)大型數(shù)據(jù)集進(jìn)行訓(xùn)練測(cè)試，網(wǎng)絡(luò)結(jié)構(gòu)及超參數(shù)定義如表1所示.網(wǎng)絡(luò)輸入為侵入深度、各向異性系數(shù)、侵入帶和原狀地層水平電阻率、圍巖層電阻率和目的層厚度所轉(zhuǎn)換為的二值圖像，即Input=BImage6×32×1, 輸出為陣列側(cè)向測(cè)井歸一化之后的視電阻率，即Output=RLA1×25.網(wǎng)格架構(gòu)搭建平臺(tái)為：Pycharm2020.2.3，Python3.7和Tensorflow2.0(GPU版本).為防止訓(xùn)練模型過(guò)擬合以及提升訓(xùn)練精度，引入Dropout和Padding機(jī)制，其中Dropout率為0.5，Padding選擇‘same’模式，即每層輸入輸出維度一致，進(jìn)而使得在每層輸入不丟失信息的前提下將卷積層延續(xù)，達(dá)到深化卷積層、提升訓(xùn)練精度的目的.除此之外，由于匯聚層將卷積層提取的多個(gè)特征融合，而對(duì)于本文回歸問(wèn)題，每個(gè)模型特征都是唯一的，因此，為保持模型特征完整性，不使用匯聚層.同時(shí)反向傳播算法采用自適應(yīng)矩估計(jì)算法(Adam).訓(xùn)練過(guò)程中，隨機(jī)選取95%的數(shù)據(jù)集為訓(xùn)練集，剩余5%的數(shù)據(jù)為測(cè)試集.

表1 2D-CNN超參數(shù)定義Table 1 The definition of hyperparameters in 2D-CNN

為測(cè)試2D-CNN網(wǎng)格訓(xùn)練精度，本文將其訓(xùn)練結(jié)果與BPNN對(duì)比分析，其中BPNN網(wǎng)絡(luò)架構(gòu)與其類似，具體為：九層神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)，激活函數(shù)為ReLU，隱含層神經(jīng)元個(gè)數(shù)分別為32，64，128，256，516，516，516，Dropout率為0.5，損失函數(shù)與反向傳播訓(xùn)練算法與2D-CNN一致.圖8為其訓(xùn)練精度與損失函數(shù)下降量的對(duì)比結(jié)果，綠色和黑色散點(diǎn)線分別代表2D-CNN和BPNN訓(xùn)練結(jié)果，可以看出兩者收斂速度大致相同，但訓(xùn)練精度卻相差較大，2D-CNN在100次迭代后訓(xùn)練精度即可達(dá)到99%左右，而B(niǎo)PNN為85%.除此之外，BPNN容易陷入局部極小值，精度無(wú)法隨著迭代的進(jìn)行而提升，相較而言，2D-CNN不僅可以快速收斂，而且訓(xùn)練精度隨著迭代的進(jìn)行逐漸接近100%.測(cè)試集與之規(guī)律類似，不再贅述.另外，通過(guò)計(jì)算2D-CNN和BPNN在測(cè)試集上的泛化誤差，計(jì)算公式如式(17).我們發(fā)現(xiàn)，兩者泛化誤差分別為1.28%和16.37%，從而也佐證了2D-CNN算法的高效性.

圖8 神經(jīng)網(wǎng)絡(luò)訓(xùn)練結(jié)果 (a) BPNN和2D-CNN預(yù)測(cè)精度； (b) BPNN和2D-CNN損失函數(shù).Fig.8 Neural network training results (a) Prediction accuracy of BPNN and 2D-CNN； (b) Loss function of BPNN and 2D-CNN.

(17)

式中，ntest為測(cè)試集樣本數(shù)，此處為211200；RLAij_2D-CNN為2D-CNN計(jì)算結(jié)果，RLAij_3D-FEM為3D-FEM計(jì)算結(jié)果.

除此之外，為驗(yàn)證訓(xùn)練模型的精度，我們隨機(jī)選取20個(gè)樣本點(diǎn)，分別測(cè)試中間層原狀地層水平電阻率與各向異性系數(shù)對(duì)模型預(yù)測(cè)精度的影響，如圖9所示，實(shí)線代表3D-FEM計(jì)算結(jié)果，散點(diǎn)表示2D-CNN預(yù)測(cè)結(jié)果.樣本點(diǎn)模型皆為三層模型：Rs為2 Ωm，H為4 m，Di為0.5 m，Rxoh為20 Ωm.其中圖9aλ為1.5，Rth變化范圍為1～500 Ωm；圖9b中Rth為100 Ωm，λ變化范圍為1～3.從對(duì)比結(jié)果可以看出，兩者吻合度極高，最大相對(duì)誤差分別為0.81%和0.24%.因此，基于深層卷積網(wǎng)絡(luò)的快速正演模型滿足精度需求，為斜井各向異性地層的快速反演提供有力支撐.

圖9 2D-CNN精度測(cè)試 (a) Rth對(duì)預(yù)測(cè)精度的影響； (b) λ對(duì)預(yù)測(cè)精度的影響.Fig.9 Prediction accuracy testing of 2D-CNN (a) The effect of Rth on the prediction accuracy； (b) The effect of λ on the prediction accuracy.

3 混合MPGA-LM算法及快速反演流程

3.1 代價(jià)函數(shù)定義

大斜度井/水平井陣列側(cè)向測(cè)井是典型的非線性反演問(wèn)題，確定性算法L-M算法由于結(jié)合高斯牛頓算法和最速下降法的特性，能夠彌補(bǔ)全局尋優(yōu)能力不強(qiáng)或收斂性差的缺陷，反演過(guò)程中采用的代價(jià)函數(shù)如下：

(18)

式中，x為待反演參數(shù)矢量，x=(Di,Rxoh,Rth,λ,Rs)T，上標(biāo)T表示矩陣的轉(zhuǎn)置，S(x)為陣列側(cè)向儀器正演響應(yīng)，dobs為實(shí)際測(cè)量數(shù)據(jù)，‖·‖2表征L2范數(shù)，xp為模型已知參考矢量，通常是上一步迭代參數(shù)矢量，ξ為正則化參數(shù)，最后一項(xiàng)的作用是壓制測(cè)量噪聲，同時(shí)減小非線性反演出現(xiàn)的病態(tài)問(wèn)題.反演的目的是得到一個(gè)x*使得代價(jià)函數(shù)最小，即x*=argmin{C(x)}.

進(jìn)一步，經(jīng)二階泰勒近似，可得第k步迭代步長(zhǎng)xk為

xk=-JT(xk){S(xk)-dobs)}[JT(xk)J(xk)+ξkI]-1,

(19)

式中，I為單位矩陣，J為S在xk處的雅克比矩陣，ξk為第k次迭代的正則化參數(shù)大小.

為克服傳統(tǒng)差分方法計(jì)算雅可比矩陣耗時(shí)的缺陷，本文采用Broyden近似方法進(jìn)行雅可比矩陣的迭代更新(Broyden, 2000).

·(xk)T,

(20)

式中，J(xk+1)為第k+1步的雅可比矩陣.

3.2 基于MPGA的初值優(yōu)選策略

由于L-M算法對(duì)初值依賴性較強(qiáng)，若初值選取不適，算法易困于局部極小值，反演精度低.為此，本文基于卷積神經(jīng)網(wǎng)絡(luò)快速正演模型，同時(shí)為克服標(biāo)準(zhǔn)遺傳算法(GA)早熟早收斂的缺陷，利用多種群遺傳算法(MPGA)優(yōu)化初值選取，進(jìn)而提升L-M算法全局尋優(yōu)能力.

如圖10所示，多種群遺傳算法的核心在于將多個(gè)種群進(jìn)行初始化，即“二進(jìn)制”編碼，并同時(shí)進(jìn)行選擇、交叉和變異操作.種群間通過(guò)移民算子進(jìn)行協(xié)同進(jìn)化，最后通過(guò)適應(yīng)度函數(shù)分析，借助人工選擇算子選取種群最優(yōu)個(gè)體組成精華種群.本文定義適應(yīng)度函數(shù)如下：

圖10 多種群遺傳算法流程Fig.10 The framework of multi-population genetic algorithm

(21)

式中，p= 1,2,…,P，P為種群數(shù)量，n=1,2,…,N，N為每個(gè)種群個(gè)體數(shù)量，m為每個(gè)測(cè)量點(diǎn)數(shù)據(jù)量，此處為25.

除此之外，本文所采用的多種群遺傳算法超參數(shù)選擇如下：種群數(shù)目為5，個(gè)體數(shù)目為30，最大遺傳代數(shù)為50，變量的二進(jìn)制位數(shù)為32，交叉概率在[0.7，0.9]內(nèi)隨機(jī)產(chǎn)生，變異概率在[0.001，0.05]內(nèi)隨機(jī)產(chǎn)生.

為測(cè)試初值選取精度，以30°三層地層模型為例，兩側(cè)為各向同性圍巖，中間層為有侵各向異性地層.表2為測(cè)試模型參數(shù)分布以及多種群遺傳算法初值選取結(jié)果，Mode l1-Mode l5層厚分別為1～5 m，間隔為1 m.可以看出，受層厚影響，當(dāng)層厚小于3 m時(shí)，圍巖電阻率可準(zhǔn)確選取，平均相對(duì)誤差為4.4%；當(dāng)層厚大于3 m時(shí)，陣列側(cè)向測(cè)井響應(yīng)基本不反映圍巖信息，選取誤差急劇增大，而此時(shí)其他參數(shù)選取較為精確，整體而言，侵入深度平均相對(duì)誤差為18.6%，各向異性系數(shù)為17.8%，侵入帶和原狀地層水平電阻率分別為18.7%和6.8%，選取精度滿足L-M算法需求.

表2 初值選取結(jié)果Table 2 The selection result of initial value

3.3 混合MPGA-LM快速反演流程及有效性測(cè)試

實(shí)際資料處理過(guò)程中，針對(duì)帶有深度點(diǎn)的大型測(cè)井響應(yīng)數(shù)據(jù)集，我們將其由上到下分為若干子區(qū)域，每個(gè)子區(qū)域?yàn)槿龑幽Ｐ?，中間層為有限厚，地層厚度由層界面位置決定，而層界面位置可通過(guò)視電阻率曲線拐點(diǎn)或伽馬、自然電位等巖性曲線確定，上下地層視為無(wú)限厚地層進(jìn)行處理.圖11為混合MPGA-LM反演流程圖，由上節(jié)可知，經(jīng)MPGA優(yōu)選后的參數(shù)分布于反演代價(jià)函數(shù)真實(shí)解附近，故可將其作為L(zhǎng)-M算法的初值，進(jìn)而提升反演效率；另一方面，由于當(dāng)層厚大于3 m時(shí)候，圍巖對(duì)陣列側(cè)向響應(yīng)近乎無(wú)影響，此時(shí)反演參數(shù)可降為4個(gè)，即Di、Rxoh、Rth和λ.同時(shí)，為提高訓(xùn)練精度，反演可循環(huán)N次(本文N為2)，即處理完整個(gè)井段資料后，可進(jìn)行二次處理.

圖11 混合MPGA-LM反演流程Fig.11 The workflow of hybrid MPGA-LM scheme

另外，為測(cè)試混合MPGA-LM快速反演算法的有效性，本文以三層模型為例，分別采用簡(jiǎn)單初始模型與MPGA提供的初始模型進(jìn)行反演對(duì)比，反演算法皆為L(zhǎng)M算法.地層傾角選擇30°，模型參數(shù)分布如表3所示.其中，在簡(jiǎn)單初始模型中分別將深、淺測(cè)井響應(yīng)RLA1和RLA5作為Rxoh和Rth的初始值，Di和λ的初始值分別為0.5 m和1.5，表3為兩者方法反演結(jié)果及相對(duì)誤差，可以看出，較之簡(jiǎn)單初始模型的LM算法，MPGA-LM可以更好的將地層剖面重構(gòu)，反演精度更高.以中間地層反演結(jié)果為例，基于簡(jiǎn)單初始模型的LM算法反演后Di、λ、Rxoh和Rth的相對(duì)誤差(RE)分別為36.67%、33%、84.10%和23.08%；MPGA-LM反演后Di、λ、Rxoh和Rth的相對(duì)誤差分別為3.33%、0.50%、6.40%和2.86%，精度可以提高數(shù)倍至數(shù)十倍，由此本文提出的混合MPGA-LM快速反演算法的有效性得到了驗(yàn)證，為下文數(shù)值模擬實(shí)例的分析提供了保障.

表3 LM與MPGA-LM算法反演結(jié)果對(duì)比Table 3 The inverted results comparison between LM and MPGA-LM algorithm

4 數(shù)值模擬實(shí)例

本節(jié)通過(guò)建立21層俄克拉荷馬模型，進(jìn)一步驗(yàn)證上述反演方案的有效性.模型參數(shù)如下：井徑為8 inch，泥漿電阻率為0.1 Ωm，地層相對(duì)傾角分別

為0°、30°和60°，侵入深度范圍0～1.0 m，侵入帶水平電阻率分布范圍為0.2～30 Ωm，原狀地層水平電阻率范圍為0.2～100 Ωm，各向異性系數(shù)范圍為1～2.5 Ωm.圖12為該模型參數(shù)分布示意圖，縱坐標(biāo)為真實(shí)垂直深度(TVD).圖13為地層相對(duì)傾角分別為0°、30°和60°時(shí)陣列側(cè)向測(cè)井響應(yīng).本文將地層界面位置作為已知條件，僅反演侵入深度、侵入帶電阻率以及原狀地層電阻率.

圖12 俄克拉荷馬模型參數(shù)分布 (a) 侵入深度； (b) 侵入帶和原狀地層水平電阻率； (c) 各向異性系數(shù).Fig.12 The Oklahoma model parameters (a) Di； (b) Rxoh and Rth；(c) λ.

圖13 陣列側(cè)向測(cè)井響應(yīng)結(jié)果 (a) 0°； (b) 30°； (c) 60°.Fig.13 The responses of array laterolog

從測(cè)井響應(yīng)結(jié)果可以看出，泥漿濾液侵入對(duì)其影響嚴(yán)重，如模型第八層，該層為各向同性地層，厚度為5.2 m，侵入深度為0.3 m，地層真實(shí)電阻率為70 Ωm，而視電阻率值在0°、30°和60°時(shí)與地層電阻率的最大偏離度，即max(RLA/Rth)，分別為32.16%、31.81%和31.03%，遠(yuǎn)遠(yuǎn)偏離地層真實(shí)電阻率；除此之外，由于受地層各向異性的影響，導(dǎo)致響應(yīng)結(jié)果介于水平和垂直電阻率之間，如模型第16層，該層無(wú)侵入，地層水平和垂直電阻率分別為5 Ωm和20 Ωm，測(cè)井響應(yīng)在0°、30°和60°時(shí)的視電阻率大小分別在12 Ωm、13 Ωm 和15 Ωm附近，可以看出，60°地層對(duì)各向異性更為敏感，由此驗(yàn)證了2.3部分的結(jié)論.由上述可知，單單從測(cè)井曲線無(wú)法提取地層真實(shí)電阻率，因此需借助反演手段對(duì)地層真實(shí)電阻率剖面進(jìn)行重構(gòu).本文基于運(yùn)行內(nèi)存為8 GB 的英特爾i7-9700 處理器對(duì)60 m測(cè)井?dāng)?shù)據(jù)進(jìn)行反演處理，數(shù)據(jù)采樣間隔為0.1 m，循環(huán)反演次數(shù)為2次，反演總耗時(shí)約6.8 min.圖14是基于本文提出的反演策略對(duì)測(cè)井資料處理后的結(jié)果.

從圖14反演結(jié)果及圖15反演誤差曲線可以看出，就反演精度而言，地層水平電阻率>侵入深度>各向異性系數(shù)，三者的相對(duì)誤差范圍分別為0～8.38%、0.11%～6.98%和0～14.35%，該現(xiàn)象可由1.3節(jié)結(jié)論解釋，即測(cè)井響應(yīng)對(duì)侵入深度的敏感性大于地層各向異性；另一方面，地層傾角對(duì)侵入深度和地層水平電阻率的反演精度影響較小，而對(duì)地層各向異性影響較大，可以看出，三個(gè)角度下各向異性的反演相對(duì)誤差最大值分別為14.35%、7.79%和3.85%，該現(xiàn)象亦可由1.3節(jié)結(jié)論解釋，即高角度條件下陣列側(cè)向測(cè)井響應(yīng)對(duì)各向異性的敏感性較強(qiáng).

圖14 混合MPGA-LM反演結(jié)果 (a) 侵入深度； (b) 侵入帶和原狀地層水平電阻率； (c) 各向異性系數(shù).Fig.14 The inverted results using hybrid MPGA-LM algorithm (a) Di； (b) Rxoh and Rth； (c) λ.

圖15 反演結(jié)果相對(duì)誤差 (a) 侵入深度反演誤差； (b) 侵入帶和原狀地層水平電阻率反演誤差； (c) 各向異性系數(shù)反演誤差.Fig.15 The relative error of inverted results (a) The relative error of inverted Di； (b) The relative error of inverted Rxoh and Rth； (c) The relative error of inverted λ.

另外，由于在反演過(guò)程中調(diào)用的正演模型是2D-CNN快速計(jì)算模型，而非3D-FEM算法，而2D-CNN模型缺乏對(duì)各向異性圍巖層的考慮.因此，為進(jìn)一步分析反演誤差產(chǎn)生原因，我們將2D-CNN與3D-CNN的計(jì)算結(jié)果進(jìn)行誤差分析，平均相對(duì)誤差(MRE)計(jì)算公式如式(22)所示，具體結(jié)果如表4所示.

表4 2D-CNN與3D-FEM計(jì)算結(jié)果平均相對(duì)誤差Table 4 The mean relative error between the responses calculated by 2D-CNN and 3D-FEM

(22)

式中，RLAi_2D-CNN為2D-CNN計(jì)算結(jié)果，RLAi_3D-FEM為3D-FEM計(jì)算結(jié)果.

可以發(fā)現(xiàn)，整體而言，兩者平均相對(duì)誤差小于10%，對(duì)于大多數(shù)層而言，相對(duì)誤差小于3%，滿足精度要求.然而在薄層處，由于受到圍巖層厚、侵入和各向異性的影響，平均相對(duì)誤差普遍大于5%，如第2層、第10層和第18層等.該現(xiàn)象正是造成反演結(jié)果精度較低的原因之一.為此，在后續(xù)研究中，我們將繼續(xù)完善訓(xùn)練模型，提升神經(jīng)網(wǎng)絡(luò)算法對(duì)圍巖層特征的拾取能力.

值得注意的是，本文所采用的反演模型均為侵入模型，由于多解性問(wèn)題的存在，會(huì)使得無(wú)侵地層的侵入深度反演結(jié)果不為0，但侵入帶電阻率和原狀地層電阻率的反演結(jié)果一致性較高，如模型第五層和第九層.整體而言，該方案可有效提升算法克服反演多解性的能力，進(jìn)而提高傾斜各向異性地層陣列側(cè)向測(cè)井反演精度，給陣列側(cè)向測(cè)井實(shí)際資料處理提供一種新思路.

5 結(jié)論

針對(duì)傾斜各向異性地層，本文基于3D-FEM，分析了侵入深度、各向異性以及地層傾角等因素對(duì)陣列側(cè)向測(cè)井響應(yīng)的影響；接著，引入2D-CNN網(wǎng)絡(luò)架構(gòu)，針對(duì)侵入深度小于1.0 m，各向異性系數(shù)小于4，層厚小于5 m，圍巖和地層水平電阻率位于0.1～500 Ωm的地層，通過(guò)建立大型數(shù)據(jù)集以及將模型可視化，使得2D-CNN訓(xùn)練精度和計(jì)算速度能夠滿足正反演需求，訓(xùn)練精度達(dá)到99%左右，計(jì)算一個(gè)測(cè)井點(diǎn)僅需0.36 ms；在此基礎(chǔ)上，通過(guò)引入全局尋優(yōu)能力較強(qiáng)的MPGA算法，可為L(zhǎng)-M反演算法提供較為精確的初值，進(jìn)而基于區(qū)域分解技術(shù)使得斜井各向異性地層電阻率反演精度達(dá)到90%以上.

不可否認(rèn)的是，在更為復(fù)雜的三維地層條件下，例如碳酸鹽巖地層，存在裂縫和溶洞發(fā)育、地層各向異性強(qiáng)以及電阻率大等特征，本文提供的數(shù)據(jù)庫(kù)及網(wǎng)絡(luò)架構(gòu)將不再適用，讀者僅需將其進(jìn)行擴(kuò)充與完善即可.另外，復(fù)雜三維地層的測(cè)井響應(yīng)會(huì)使得反演多解性更為強(qiáng)烈，為此，在今后研究中，應(yīng)進(jìn)一步豐富和完善神經(jīng)網(wǎng)絡(luò)架構(gòu)以及反演方案，以期為油氣田開(kāi)發(fā)提供更為便捷可靠的處理方案.