基于改進(jìn)雙向峰-谷搜索算法的煤矸模型電磁波正演模擬

史翔予， 司壘， 王忠賓， 魏東， 顧進(jìn)恒

（中國礦業(yè)大學(xué) 機電工程學(xué)院，江蘇 徐州 221116）

0 引言

2022年，全國煤炭產(chǎn)量約45.6億t，同比增長10.5%，創(chuàng)歷史新高，煤炭消費量占全國能源消費總量的56.2%，這體現(xiàn)了煤炭在我國能源結(jié)構(gòu)中的主體地位。厚煤層儲量和產(chǎn)量在總量中的占比均為44%左右，因此厚煤層的高效開采對于我國煤炭開采行業(yè)具有重要意義。綜放開采技術(shù)是我國厚煤層開采的主要方法。我國綜采技術(shù)起步早，發(fā)展快，目前機采已經(jīng)可以實現(xiàn)自動化控制，但是由于放煤過程中煤矸含量識別技術(shù)尚未成熟，放煤部分仍只能采用人工操作，極易造成欠放、過放等情況，大大降低煤炭采出率[1-2]。因此，實現(xiàn)放頂煤過程煤矸含量自動識別是綜采自動化的重要目標(biāo)。

王家臣等[3-4]利用圖像識別技術(shù)將混矸率劃分為投影面積混矸率、表面體積混矸率、內(nèi)部體積混矸率，實現(xiàn)投影面積混矸率向表面體積混矸率的過渡，但井下環(huán)境低照度、高粉塵仍是制約圖像識別技術(shù)實際應(yīng)用的重要因素。竇希杰等[5-6]采用經(jīng)驗?zāi)B(tài)分解方法將綜放工作面煤和矸石沖擊液壓支架尾梁產(chǎn)生的振動信號分解成多個本征模態(tài)函數(shù)分量，采用多種機器學(xué)習(xí)方法訓(xùn)練識別模型，其中隨機森林模型準(zhǔn)確率最高。馬瑞等[7]、宋慶軍等[8]依據(jù)煤和矸石與刮板輸送機發(fā)生碰撞產(chǎn)生的聲波頻率的差異進(jìn)行煤矸識別，但由于工作面存在大量噪聲，影響結(jié)果準(zhǔn)確性。李瑞等[9]、丁震等[10]利用紅外光譜技術(shù)結(jié)合黑白校正、SG卷積平滑等預(yù)處理方法建立了煤矸分類模型，對不同礦區(qū)的煤與矸石都有良好的適用性，但是實時性較差，無法滿足工程需要。張寧波等[11]、王增才[12]根據(jù)γ射線穿透煤矸混合物時發(fā)生的衰減規(guī)律，建立了煤矸混合體灰分與含矸量間的量化關(guān)系，但是不同的巖石包含放射性含量也不同，影響識別精度。

電磁波波速快，穿透性強，可以避免受到粉塵及光纖干擾的影響；電磁波攜帶的信息豐富，在不同電性參數(shù)的介質(zhì)中所體現(xiàn)的傳播特性也不同。放頂煤過程產(chǎn)生的煤矸混合物是由煤、矸石和空氣形成的三相介質(zhì)，各相介質(zhì)的電性參數(shù)不同、空間位置分布隨機且混合形式復(fù)雜。在不同組分的混合三相介質(zhì)中，電磁波的傳播特性也不盡相同。煤和矸石的相對介電常數(shù)差異明顯，通過研究不同含矸率煤矸混合物的電性參數(shù)，可為放頂煤工作面含矸率自動識別提供新的思路和方法。為了探究不同含矸率煤矸混合物的電性差異，本文提出了一種基于分治策略的雙向峰-谷搜索算法，建立了煤矸多相離散隨機介質(zhì)模型，基于麥克斯韋方程組及其本構(gòu)關(guān)系方程，利用時域有限差分法對所建立的模型進(jìn)行電磁波正演模擬。

1 煤矸多相離散隨機介質(zhì)模型

綜放開采形成的煤矸混合物隨機性較強，傳統(tǒng)的均勻介質(zhì)模型無法精準(zhǔn)描述煤-矸石-空氣的三相隨機分布狀態(tài)，無法為電磁波探測相關(guān)研究提供滿足需求的樣本。因此，需要建立能夠準(zhǔn)確描述煤矸混合物空間分布和電性參數(shù)的介質(zhì)模型，在此基礎(chǔ)上對其進(jìn)行電磁波正演仿真，實現(xiàn)電性參數(shù)差異分析。

1.1 隨機介質(zhì)理論

基于隨機介質(zhì)理論的隨機介質(zhì)模型最早出現(xiàn)于地質(zhì)學(xué)與大地構(gòu)造物理學(xué)的相關(guān)研究中。地質(zhì)體的相組成非常復(fù)雜，為非均介質(zhì)，其中包含大量孔隙、斷巖、物質(zhì)夾雜等。從相結(jié)構(gòu)來分析，地質(zhì)體的相結(jié)構(gòu)與放頂煤形成的煤-矸石-空氣多相結(jié)構(gòu)較為相似。隨機介質(zhì)可以分解為大、小2種尺度的非均勻性部分：大尺度非均性用于描述介質(zhì)的背景特征，如層狀連續(xù)介質(zhì)模型；小尺度非均性則用于描述背景特征上的隨機擾動，通常用一個均值為零的二階平穩(wěn)隨機過程來表示[13]。

對于煤矸混合物隨機介質(zhì)模型，其相對介電常數(shù)為

式中：εm(x,y)為隨機介質(zhì)的大尺度非均勻性部分，與坐標(biāo)(x,y)有一定關(guān)系；δm(x,y)為隨機介質(zhì)的小尺度非均勻性部分。

δm(x,y)在二階平穩(wěn)假設(shè)下可分解為

式中：δp(x,y)為相對介電常數(shù)變化的標(biāo)準(zhǔn)差；f(x,y)為均值為零、標(biāo)準(zhǔn)差為1的二維隨機序列，其空間分布服從二維空間自相關(guān)函數(shù)。

為了靈活準(zhǔn)確地描述實際介質(zhì)，二維空間自相關(guān)l函數(shù)通常采用高斯-指數(shù)混合型橢圓自相關(guān)函數(shù)[14]：

式中：i，j分別為隨機介質(zhì)的自相關(guān)長度和均勻性；r為粗糙度因子。

1.2 基于分治策略的雙向峰-谷搜索算法

隨機介質(zhì)模型是一種連續(xù)的空間分布狀態(tài)，且相對介電常數(shù)隨機變化[15]，然而綜放開采形成的煤矸混合物是由相對介電常數(shù)存在絕對差異的煤炭、矸石、空氣構(gòu)成，各組成之間存在明顯的相界面且滿足一定體積比，因此，必須采用合適的方法進(jìn)行離散性處理。

Lin Li等[16]首先提出將雙向峰-谷搜索算法用于建立鋁硅聚苯酯（AlSi-PHB）封嚴(yán)涂層的多相離散隨機介質(zhì)模型。其將隨機介質(zhì)理論與統(tǒng)計學(xué)方法引入封嚴(yán)涂層建模研究，用大尺度非均勻性描述涂層的平均特性，在該平均值上添加隨機擾動，實現(xiàn)小尺度非均勻性描述，然后采用雙向峰-谷搜索算法對隨機介質(zhì)模型進(jìn)行離散性處理，處理后的模型在空間上離散分布且各組分滿足一定體積比，能準(zhǔn)確描述封嚴(yán)涂層微觀形貌。綜放形成的煤矸混合物與AlSi-PHB封嚴(yán)涂層存在多相、非均質(zhì)、組成形貌隨機、各組分物理性質(zhì)差異大等共同特點。AlSi與PHB存在絕對密度差異，而煤和矸石存在相對介電常數(shù)的差異，因此，本文利用雙向峰-谷搜索算法建立煤-矸石-空氣三相介質(zhì)的多相離散隨機介質(zhì)模型。

相對介電常數(shù)的差異源于物質(zhì)的化學(xué)組成?？諝獾南鄬殡姵?shù)一般認(rèn)為是1。煤的主要成分是有機物，儲存電荷的能力相對較弱，屬于弱極性物質(zhì)，不同類型的煤的化合物成分有細(xì)微差別，其相對介電常數(shù)介于3和3.8之間。矸石的主要成分是無機化合物，儲存電荷能力強，屬于極性物質(zhì)，根據(jù)相關(guān)文獻(xiàn)[17-18]，取矸石相對介電常數(shù)為7。

用雙向峰-谷搜索算法建立具有指定含矸率、含煤率的多相離散隨機介質(zhì)模型的步驟如下：

1） 確定矸石、煤及空氣的相對介電常數(shù)ε1—ε3。

2） 根據(jù)介質(zhì)中相對介電常數(shù)的分布，沿最大值方向搜索到根據(jù)含矸率設(shè)置的閾值后，將此部分設(shè)置為矸石相；同時沿最小值方向搜索到根據(jù)空氣含量設(shè)置的閾值后，將此部分設(shè)置為空氣相。

3） 遍歷網(wǎng)格中所有點，判斷其屬性值是否為矸石的相對介電常數(shù)ε1或空氣的相對介電常數(shù)ε3，若是則不處理，若不是則將該點相對介電常數(shù)設(shè)置為煤炭的相對介電常數(shù)ε2。

以含矸率30%、空氣含量10%、含煤率60%為例，建立多相離散隨機介質(zhì)模型，如圖1所示，深藍(lán)色部分是空氣，淺藍(lán)色部分是煤炭，黃色部分是矸石。

圖1 算法改進(jìn)前多相離散隨機介質(zhì)模型Fig. 1 Multiphase discrete random medium model before algorithm improvement

由圖1可知，煤、矸石、空氣三相之間存在明顯的相界面。然而，各組分在模型區(qū)域內(nèi)呈現(xiàn)聚集分布，局部區(qū)域呈現(xiàn)大量空隙或大量煤炭聚集的現(xiàn)象，導(dǎo)致部分局部無法表征整體的電性參數(shù)。例如，從縱向網(wǎng)格線150處入射的電磁波，其橫向穿透區(qū)域大部分是相對介電常數(shù)較低的空氣；而從縱向網(wǎng)格線100處入射的電磁波，其橫向穿透區(qū)域的煤炭體積分?jǐn)?shù)超過80%，從而引起電磁波正演仿真的失真。針對該問題，本文引入分治策略對雙向峰-谷搜索算法進(jìn)行改進(jìn)。

分治策略應(yīng)用遞歸的思想，首先將問題劃分為子問題，然后對子問題遞歸求解，最后將子問題的解合并成原問題的解。雙向峰-谷搜索算法在設(shè)定閾值后直接對整個模型空間進(jìn)行建模，而基于分治策略改進(jìn)的算法步驟如下：

1） 將整個模型等幅劃分為100個子網(wǎng)格并編號。

2） 通過以下方法求解子網(wǎng)格順序：計算子網(wǎng)格中的含煤率p1和含矸率φ1，并確定歸一化系數(shù)（p，φ分別為設(shè)定的含煤率和含矸率），將其代入傅里葉歸一化方程中進(jìn)行迭代。

3） 每次迭代后依據(jù)更新后的歸一化系數(shù)反推出p1和φ1，并將其作為模型閾值重建子網(wǎng)格中的煤矸模型，直至歸一化系數(shù)迭代到1，此時該網(wǎng)格內(nèi)的煤矸含量與設(shè)定值相同，子網(wǎng)格求解完成。

4） 將各個子網(wǎng)格的解按順序合并，完成整個模型空間的建模。

改進(jìn)雙向峰-谷搜索算法流程如圖2所示，其中虛線部分表示算法執(zhí)行多次迭代。

圖2 改進(jìn)雙向峰-谷搜索算法流程Fig. 2 Flow of improved bidirectional peak-valley search algorithm

傅里葉歸一化方程可通過迭代令參數(shù)平滑收斂于1，通過更改自設(shè)定常數(shù)來調(diào)節(jié)收斂速度，其一般形式為

式中：Γ*為迭代更新后的歸一化系數(shù)；a0—a4，w，b1—b4為自設(shè)定常數(shù)；Γ為原歸一化系數(shù)。

同樣，以含矸率30%、含煤率60%為例，用改進(jìn)后的算法建立煤矸多相離散隨機介質(zhì)模型，如圖3所示。

圖3 算法改進(jìn)后煤矸多相離散隨機介質(zhì)模型Fig. 3 Multiphase discrete random medium model after algorithm improvement

為了比較算法改進(jìn)前后的性能，本文運用二值法，根據(jù)所建立模型各相的灰度值不同，將算法改進(jìn)前后模型的三相單獨分離出來，然后利用ImageJ軟件基于面積分?jǐn)?shù)進(jìn)行分析。算法改進(jìn)前后模型各相分布如圖4所示。

圖4 算法改進(jìn)前后模型各相分布Fig. 4 Each phase distribution of model before and after algorithm improvement

分析結(jié)果如下：改進(jìn)前模型空氣相占10.28%，煤相占59.72%，矸石相占30.00%；改進(jìn)后模型空氣相占10.16%，煤相占59.84%，矸石相占30.00%?？紤]到基于灰度值差異的二值分離法易在空氣相與煤相交叉處產(chǎn)生混雜，從而存在一定誤差，但二者誤差皆在允許范圍內(nèi)，因此算法改進(jìn)前后模型均滿足設(shè)定體積比。然而從圖4可看出，相較于算法改進(jìn)前模型，算法改進(jìn)后模型中的煤、矸石、空氣三相之間不僅存在明顯的相界面，且各相離散程度更大，不存在聚集現(xiàn)象，因此局部介質(zhì)也能體現(xiàn)整體的電性參數(shù)，能夠滿足電磁波正演的介質(zhì)模型需求。

2 電磁波正演模擬基礎(chǔ)理論

2.1 電磁波基礎(chǔ)理論

在宏觀尺度上，電磁波的產(chǎn)生與傳播都可以通過麥克斯韋方程描述，因此電磁波正演模擬的一切電磁過程都可以從麥克斯韋方程及其本構(gòu)關(guān)系方程出發(fā)進(jìn)行分析。麥克斯韋方程組描述了變化的磁場與電場的感應(yīng)關(guān)系及變化的電流與磁場的感應(yīng)關(guān)系，其微分形式為[19]

式中：?為哈密頓算子；H為磁場強度，A/m；J為電流密度，A/m2；D為電位移，C/m2；t為時間，s；E為電場強度，V/m；B為磁感應(yīng)強度，T；ρ為電荷密度，C/m3。

本構(gòu)關(guān)系方程描述了電磁場中介質(zhì)的電磁參數(shù)對電磁場的影響，其形式為

式中：ε為介質(zhì)的相對介電常數(shù)，F(xiàn)/m；μ為介質(zhì)的磁導(dǎo)率，H/m；σ為介質(zhì)的電導(dǎo)率，S/m。

2.2 時域有限差分法

時域有限差分法最早于1996年被提出，用來計算電磁波傳導(dǎo)問題。時域有限差分法是基于麥克斯韋方程的時域微分算法，其依據(jù)Yee網(wǎng)格（圖5）中各場量（電場Ex，Ey，Ez，磁場Hx，Hy，Hz）分布，將麥克斯韋方程從微分形式轉(zhuǎn)換為差分形式，然后以該網(wǎng)格單元將仿真區(qū)間離散劃分，在微小網(wǎng)格單元上采用步進(jìn)法在時間上迭代計算連續(xù)電磁場。在計算過程中通過電場與磁場相互更新，從而實現(xiàn)電磁波的迭代計算。

圖5 時域有限差分法Yee網(wǎng)格Fig. 5 Yee grid of finite difference time domain method

在實際應(yīng)用時，時域有限差分法空間及時間步長的設(shè)定會引起數(shù)值色散，從而使波形產(chǎn)生畸變和各向異性失真。然而由于計算機內(nèi)存空間不能取無限小，所以數(shù)值色散不可避免。為了減小數(shù)值色散的影響，一般采取經(jīng)驗公式進(jìn)行空間和步長選取[20]：

式中：?x,?y,?z分別為空間坐標(biāo)x，y，z的變化量；λmin為計算頻段內(nèi)電磁波在仿真空間各介質(zhì)中的最小波長。

當(dāng)?shù)胶涂臻g步長確定后，時間步長可以由庫朗（Courant）條件確定。庫朗條件要求波在網(wǎng)格中的傳播速度小于迭代所允許的最大速度，因此時間步長?t的限制條件為[21]

式中c為自由空間中的光速，m/s。

3 電磁波正演模擬

為了探究放頂煤形成的煤矸混合物含矸率對電磁波傳播的影響，最經(jīng)濟(jì)且直觀的方法是建立不同含矸率的煤矸模型進(jìn)行電磁波正演模擬，此方法的重點在于所建立的模型能夠準(zhǔn)確描述不同含矸率煤矸混合物的電性參數(shù)。最傳統(tǒng)的均勻介質(zhì)模型只適用于建立理想狀態(tài)下的介質(zhì)模型，如層狀連續(xù)介質(zhì)模型；使用未改進(jìn)的雙向峰-谷搜索算法建立的多相離散隨機介質(zhì)模型雖然能夠描述煤矸混合物中各相隨機分布的空間狀態(tài)，但是電性參數(shù)呈現(xiàn)局部突變性，導(dǎo)致電磁波正演仿真結(jié)果出現(xiàn)隨機性，無法展現(xiàn)真實規(guī)律；而基于分治策略對雙向峰-谷搜索算法進(jìn)行改進(jìn)后，多相離散隨機介質(zhì)模型不僅能正確描述煤矸混合物各相隨機分布的空間狀態(tài)，還能準(zhǔn)確刻畫其電性參數(shù)，使用此模型進(jìn)行電磁波正演仿真能夠展現(xiàn)電磁波在不同含矸率煤矸混合物中傳播的真實規(guī)律。

3.1 模擬環(huán)境

電磁波正演模擬如圖6所示。激勵平面發(fā)出一定頻率的電磁波信號，接收平面接收電磁波的透射波信號，模擬區(qū)域網(wǎng)格空間為600×600，網(wǎng)格步長為0.4。由于仿真空間有限，電磁波會在仿真空間內(nèi)發(fā)生多次反射，產(chǎn)生的反射波會擾亂接收界面的數(shù)據(jù)接收，導(dǎo)致仿真結(jié)果失真。為了使仿真條件更加貼合實際工況，模擬電磁波在遠(yuǎn)大于仿真區(qū)間的區(qū)域內(nèi)進(jìn)行傳播，本文在仿真區(qū)間的邊界處設(shè)置完美匹配層（Perfectly Matched Layer，PML）邊界條件，用來吸收邊界處的反射電磁波。

圖6 電磁波正演模擬Fig. 6 Electromagnetic wave forward simulation

在電磁波探測領(lǐng)域，激勵源常用雷克子波（Ricker），其具有延續(xù)時間短、收斂迅速、波形簡單等優(yōu)點。激勵源的激勵頻率一般取12 GHz以下，激勵頻率過高會導(dǎo)致電磁波從導(dǎo)體的外圍通過，產(chǎn)生“趨膚效應(yīng)”。為確定正演模擬的激勵頻率，本文以圖3所示30%含矸率的多相離散隨機介質(zhì)模型作為試驗介質(zhì)，分別將電磁波激勵頻率設(shè)置為2，4，6，8，10，12 GHz，從接收平面獲取的信號時域圖像如圖7所示。

圖7 不同激勵頻率下接收平面獲取的信號時域圖像Fig. 7 Time domain images of signals obtained from the receiving plane under different excitation frequencies

從圖7可看出：① 不同激勵頻率對信號透射時間影響不大，都在3.2 ns左右，這是因為煤矸混合物介質(zhì)是非色散介質(zhì)，其相對介電常數(shù)與電磁波頻率無關(guān)，而介質(zhì)對電磁波的折射率由相對介電常數(shù)決定，根據(jù)電磁波在介質(zhì)中的傳播速度計算公式，無論激勵頻率大小，電磁波在相同組分含量的煤矸混合物介質(zhì)中的傳播速度均相同。② 不同激勵頻率的信號幅值不同，其中激勵頻率為2 GHz時，信號波形出現(xiàn)了擾動，原因是激勵頻率過低時信號不穩(wěn)定，且低激勵頻率在實際工況中抗干擾性差。綜上所述，采用頻率為6～12 GHz的激勵信號可以在保證信號穩(wěn)定的前提下，信號幅值處在易于觀測的區(qū)間，且避免產(chǎn)生“趨膚效應(yīng)”。

3.2 模擬結(jié)果

用改進(jìn)算法建立6種不同含矸率（0，10%，20%，30%，40%，50%）的煤矸多相離散隨機介質(zhì)模型，如圖8所示，將其作為樣本模型，以頻率8 GHz的雷克子波作為激勵源，開展電磁波正演模擬。

圖8 不同含矸率的煤矸多相離散隨機介質(zhì)模型Fig. 8 Multiphase discrete random medium model of coal gangue with different gangue content

使用時域有限差分法對煤矸多相離散隨機介質(zhì)模型進(jìn)行電磁波正演模擬，通過解算麥克斯韋方程，得到透射波時域波形，如圖9所示。

圖9 透射波時域波形Fig. 9 Time domain waveforms of transmitted wave

從圖9可以看出：透過含矸率為0的煤矸模型的電磁波信號最先達(dá)到峰值，且信號幅值最高；透過含矸率為50%的煤矸模型的電磁波信號最后達(dá)到峰值，且信號幅值最低；隨著含矸率升高，透射電磁波的信號幅值依次降低，并且透射時間依次延長。由于矸石的相對介電常數(shù)遠(yuǎn)高于煤和空氣，所以含矸率升高使煤矸混合物整體的等效介電常數(shù)升高。這在一方面會引起介質(zhì)整體的電導(dǎo)率增加，而電磁波在介質(zhì)中的損耗與電導(dǎo)率呈正比，故電磁波在含矸率越高的模型中損耗越大，信號幅值越低；另一方面會引起介質(zhì)的折射率增加，電磁波的傳播速度與折射率成反比，故電磁波在含矸率越高的模型中透射時間越長。含矸率為40%的模型與含矸率為50%的模型的峰值時間和信號幅值差異不大，這是由于含矸率達(dá)到一定值后，模型整體的相對介電常數(shù)趨于巖石的相對介電常數(shù)。工程上含矸率通?？刂圃?0%～30%，而從圖9可看出，在此區(qū)間電磁波的傳播特性差異明顯，因此本方法具備實際應(yīng)用條件。

4 結(jié)論

1） 放頂煤形成的煤矸混合物不適合使用傳統(tǒng)的連續(xù)隨機介質(zhì)描述，采用分治策略改進(jìn)雙向峰-谷搜索算法，并基于改進(jìn)算法建立煤矸多相離散隨機介質(zhì)模型，更加貼合煤矸混合物的實際特點，能夠在電磁波正演中體現(xiàn)不同含矸率在電性參數(shù)上的差別。

2） 激勵信號的頻率會影響透射波的幅值：在12 GHz范圍內(nèi)，激勵信號頻率越高，透射波幅值越大；頻率過低會降低信號的魯棒性，激勵頻率應(yīng)高于2 GHz。

3） 煤矸混合物的含矸率與介質(zhì)整體的等效介電常數(shù)呈正相關(guān)。含矸率越高，電磁波信號的傳播損耗越多，接收平面接收到的信號幅值越小，電磁波信號穿透介質(zhì)所用的時間越長。不同含矸率之間呈現(xiàn)明顯的差異性，可以用作綜采放頂煤含矸率識別的依據(jù)。

4） 電磁波傳播速度極快，可在納秒級別完成信號收發(fā)，在實際工況下，短時間內(nèi)便可獲得大量數(shù)據(jù)集。這對于將來使用機器學(xué)習(xí)等人工智能技術(shù)提高煤矸識別速度和準(zhǔn)確度有重要意義，也說明將電磁波探測技術(shù)應(yīng)用于綜放工作面智能放煤過程中含矸率自動識別領(lǐng)域有較好的發(fā)展前景。