付宇,艾寒冰,3,姚振岸,梅竹虛,蘇可嘉,4
(1.江西省防震減災(zāi)與工程地質(zhì)災(zāi)害探測工程研究中心,江西 南昌 330013;2.東華理工大學(xué) 地球物理與測控技術(shù)學(xué)院,江西 南昌 330013;3.中國地質(zhì)大學(xué)(武漢) 地球物理與空間信息學(xué)院,湖北 武漢 430074;4.核工業(yè)二七0研究所,江西 南昌 330200)
隨著工程建設(shè)場地復(fù)雜程度以及施工難度的提高,獲得精確的地層信息對于工程的設(shè)計(jì)和施工變得尤為重要[1]。瑞利波勘探具有淺層分辨率高、抗干擾能力強(qiáng)等優(yōu)點(diǎn),在工程勘探中應(yīng)用廣泛[2]。頻散曲線反演是瑞利波勘探的關(guān)鍵環(huán)節(jié),其反演效果直接決定了獲取地層信息的精度。因此,準(zhǔn)確地進(jìn)行頻散曲線反演極其重要。
頻散曲線反演隸屬地球物理反演問題范疇,具有高度非線性、多參數(shù)、多極值的特點(diǎn)[3-4]。同時(shí)由于地震數(shù)據(jù)存在噪聲、反演多解性等因素的影響,想要精確且高效地反演頻散曲線變得非常困難。頻散曲線的反演算法可分為局部優(yōu)化算法和全局優(yōu)化算法。局部優(yōu)化算法有最小二乘法[5-9]、L-M算法[10]、Occam算法[11-12]等。早期,局部優(yōu)化算法因計(jì)算速度快而被廣泛應(yīng)用。但此類方法過度依賴于初始模型的建立,初始模型選取和真實(shí)情況相近才能獲得較好的反演效果。另外,算法易陷入局部極值,并且涉及到的偏導(dǎo)數(shù)計(jì)算以及反演結(jié)果的好壞都受到雅可比矩陣精度的影響[13-14]。這些使局部優(yōu)化算法在瑞利波反演領(lǐng)域的發(fā)展受到限制。另一種全局優(yōu)化算法因可以很好地避免局部優(yōu)化算法在反演時(shí)遇到的問題,受到學(xué)者們的廣泛關(guān)注。全局優(yōu)化算法中常用的有遺傳算法[15-17]、模擬退火算法[18-22]、粒子群算法等[23-24]。這些算法和局部優(yōu)化算法相比雖具有不依賴于初始模型、不易陷入局部最優(yōu)解等優(yōu)勢,但存在運(yùn)算時(shí)間長、容易早熟收斂等缺點(diǎn)。
鑒于此,本文將一種新的正余弦算法(sine cosine algorithm,SCA)用于瑞利波頻散曲線的反演研究。為了評價(jià)SCA對瑞利波頻散曲線反演的可行性、有效性、穩(wěn)定性和實(shí)用性,本文首先利用SCA對4個(gè)從簡單到復(fù)雜的地質(zhì)模型產(chǎn)生的含噪聲與不含噪聲的理論頻散曲線進(jìn)行反演,并分析與評價(jià)SCA的可行性,有效性和穩(wěn)定性。接著,通過進(jìn)一步對比SCA與PSO的反演性能以說明SCA相對于經(jīng)典的粒子群算法是否具有更穩(wěn)定、收斂速度更快,反演精度更高的優(yōu)越性。最后將SCA用于反演來自冰島Arnarb?lidi和美國懷俄明地區(qū)的實(shí)測數(shù)據(jù)以檢驗(yàn)SCA對瑞利波頻散曲線反演的實(shí)用性。
SCA是澳洲學(xué)者Seyedali Mirjalili在基于正余弦函數(shù)組合,應(yīng)用多個(gè)隨機(jī)參數(shù)和自適應(yīng)變量調(diào)整尋優(yōu)過程的基礎(chǔ)上于2016年提出的一種新型群智能優(yōu)化算法[25]。SCA與其他群智能優(yōu)化算法相同,可以根據(jù)先驗(yàn)信息給出每個(gè)解的初始值。如果先驗(yàn)信息不可用,則可以隨機(jī)生成初始集合。然后這些解的位置更新方程如下:
(1)
(2)
式中:t為當(dāng)前迭代次數(shù);T為最大迭代次數(shù);a為一常量。圖1給出了SCA的算法流程。
圖1 SCA算法流程
在實(shí)際工程應(yīng)用中,基階頻散曲線能量更強(qiáng),更易觀測,采集數(shù)據(jù)中大多只含有基階信息,所以在反演時(shí)一般選取基階頻散曲線進(jìn)行反演[26]。本文考慮到出現(xiàn)高階頻散的情況,對基階、高階頻散曲線聯(lián)合反演的情況也進(jìn)行了測試。在測試中將所有地層信息都設(shè)為未知,即同時(shí)反演厚度、密度、橫波速度以及縱波速度,最后選取對頻散曲線有顯著影響的層厚和橫波速度來進(jìn)行分析與評價(jià)。由于真實(shí)的地層情況比較復(fù)雜,將模型參數(shù)的搜索范圍設(shè)置為真實(shí)模型的50%。對于算法參數(shù)設(shè)置,SCA中影響算法性能的參數(shù)主要為種群數(shù)和迭代次數(shù)。本文設(shè)置SCA的種群數(shù)為30,迭代次數(shù)為100。為避免算法的隨機(jī)性,每次理論模型測試都進(jìn)行獨(dú)立反演30次,且每次反演初始模型都不相同,最后將30次反演所得數(shù)據(jù)的均值作為反演結(jié)果進(jìn)而輸出,標(biāo)準(zhǔn)差作為衡量算法穩(wěn)定性的指標(biāo)進(jìn)而評價(jià)。
瑞利波頻散曲線反演的本質(zhì)是求解適應(yīng)度函數(shù)最小值的優(yōu)化問題。本文采用的適應(yīng)度函數(shù)是根據(jù)反演所得模型能否精確解釋觀測資料所定[27]:
(3)
為合理評價(jià)SCA算法的性能,本文使用4個(gè)工程勘察中常見的理論模型對其進(jìn)行測試。這些模型從簡單到復(fù)雜逐漸過渡,模型參數(shù)見表1。表1中模型A為兩層速度遞增模型;模型B為4層速度遞增模型;模型C為4層含一低速軟弱夾層模型;模型D為4層含一高速硬夾層模型。
表1 模型參數(shù)及反演搜索范圍
為合理評測算法,從基礎(chǔ)兩層模型開始進(jìn)行測試,隨后提升模型復(fù)雜度,在四層復(fù)雜模型中進(jìn)行進(jìn)一步評測。從圖2簡單兩層速度遞增模型反演結(jié)果可看出,在無任何先驗(yàn)信息的情況下,反演所得頻散曲線(圖2a中的點(diǎn)線)和理論模擬頻散曲線(圖2a中的實(shí)線)仍具有較高擬合度,并且反演模型和理論模型各參數(shù)誤差只有0.23%,0.60%和4.76%(表2),30次獨(dú)立反演所得參數(shù)標(biāo)準(zhǔn)差表明30次獨(dú)立反演整體偏離度較低,驗(yàn)證了SCA用于簡單地層結(jié)構(gòu)頻散曲線反演的可行性以及穩(wěn)定性。
表2 模型A、B、C、D含噪聲與不含噪聲反演結(jié)果
a—反演所得頻散曲線;b—反演的橫波速度剖面
為進(jìn)一步測試SCA在復(fù)雜地層中的反演表現(xiàn),將其應(yīng)用在3種復(fù)雜地層:模型B、 C和D中進(jìn)行反演試算。結(jié)果如圖3所示,其中a2,b2和c2分別展示了模型B、 C和D的反演結(jié)果。從反演結(jié)果可以看出,不論是3種模型中的哪一種,利用SCA反演得到的頻散曲線(圖中a1,b1和c1中的點(diǎn)線)都能與理論頻散曲線(圖中a1,b1和c1中的實(shí)線)吻合良好。3種模型反演所得參數(shù)中最大相對誤差僅為6.76%,標(biāo)準(zhǔn)差值也相對較小(表2),說明SCA在處理復(fù)雜地層模型的頻散曲線反演問題時(shí)依舊能保持良好的性能。
a1,b1,c1—分別為模型B、C和D反演所得頻散曲線;a2,b2,c2—分別為模型B、C和D反演的橫波速度剖面
在采集地震數(shù)據(jù)時(shí),噪聲是難以避免的。這些噪聲會使提取到的頻散曲線在小范圍內(nèi)隨機(jī)擾動(dòng),進(jìn)而影響反演結(jié)果。因此,算法用于實(shí)際數(shù)據(jù)處理前需檢測其抗噪能力。為測試SCA的抗噪性能,在以上4個(gè)理論模型經(jīng)正演模擬得到的頻散曲線中加入15%的隨機(jī)噪聲后進(jìn)行反演。
SCA用于含噪聲模型的反演結(jié)果如圖4所示,其中a2,b2,c2和d2分別展示了模型A、B、C、D的反演結(jié)果。從圖a2,b2,c2和d2可以看出,即便加入15%的噪聲,反演模型參數(shù)和真實(shí)模型參數(shù)依然相近,并且模型A、B、C、D利用SCA反演所得參數(shù)的相對誤差僅變?yōu)榱?.00%,6.54%,8.16%和1.36%(表2)。將其和不加噪聲模型反演結(jié)果進(jìn)行對比,可以看出反演所得參數(shù)間幾乎無偏差(表2),這表明噪聲雖然對反演結(jié)果具有一定影響,但SCA具有較強(qiáng)的抗噪能力。故SCA適用于包含較強(qiáng)噪聲數(shù)據(jù)的反演。
a1,b1,c1,d1—分別為模型A、B、C和D反演所得頻散曲線;a2,b2,c2,d2—分別為模型A、B、C和D反演的橫波速度剖面
在一些特殊地層(如軟弱夾層),可能會出現(xiàn)高階波高頻部分能量強(qiáng)于基階波的情況[28],這時(shí)若將高階頻散曲線和基階頻散曲線聯(lián)合起來反演,便可以增加反演的有效信息,進(jìn)而提升反演所得地層信息的準(zhǔn)確度,故檢驗(yàn)算法對于多模式頻散曲線的反演能力是必要的。此次測試以模型C為例,向模型C的基階頻散曲線數(shù)據(jù)中加入高階數(shù)據(jù)后進(jìn)行反演試算,反演結(jié)果見圖5和表3。從圖5a1可知,反演得到的頻散曲線(點(diǎn)線),無論基階還是高階,都與真實(shí)模型頻散曲線(實(shí)線)擬合得很好并且由圖5a2可看出,反演所得模型參數(shù)(點(diǎn)線)與真實(shí)模型參數(shù)(實(shí)線)相差也較小,說明SCA用于多模式頻散曲線反演是可行的。另外,多模式頻散曲線反演結(jié)果相比單純基階頻散曲線反演結(jié)果的精度總體有所提升,尤其是第三層層厚的相對誤差從6.28%降到了1.38%。由此可知,多模式頻散曲線反演比只用基階頻散曲線進(jìn)行反演所得模型參數(shù)準(zhǔn)確度更高,反演效果更優(yōu)異。
表3 模型C基階數(shù)據(jù)和多模式數(shù)據(jù)反演結(jié)果
a—反演所得頻散曲線;b—反演的橫波速度剖面
為證實(shí)SCA用于頻散曲線反演確實(shí)能得到較為準(zhǔn)確的地層參數(shù),將其和經(jīng)典的粒子群算法(particle swarm optimization,PSO)一起進(jìn)行反演測試,對比、分析和評價(jià)兩者的反演性能。測試在不含噪聲的模型D上進(jìn)行,各單獨(dú)反演30次。反演時(shí)兩種算法的種群數(shù)目,搜索空間和迭代次數(shù)都相同。反演結(jié)果如圖6所示,具體的反演參數(shù)如表4所示。結(jié)合圖6a2和表4可知,SCA在迭代初期就大致確定了最優(yōu)解位置,迭代到30次的時(shí)候就基本收斂且收斂精度較高為1.35 m/s;PSO的尋優(yōu)能力較差,在迭代到40次左右才收斂并且在100次左右時(shí)還有進(jìn)一步收斂的趨勢,最后收斂的精度也較低為11.57 m/s。除此之外,通過表4中SCA與PSO對于模型各參數(shù)反演的標(biāo)準(zhǔn)差可以看到,SCA算法的穩(wěn)定性遠(yuǎn)遠(yuǎn)強(qiáng)于傳統(tǒng)的粒子群算法,這得益于SCA獨(dú)特的更新解的方式。綜上可知SCA相較于PSO在瑞雷波頻散曲線反演中的表現(xiàn)具有更高的精度,更快的收斂速度以及更強(qiáng)的穩(wěn)定性。
表4 SCA和PSO反演效果對比
a—放大前的收斂曲線對比;b—放大后的收斂曲線對比
前面理論模型測試結(jié)果表明,SCA用于頻散曲線反演是可行且有效的。現(xiàn)用冰島Arnarb?lidi地區(qū)數(shù)據(jù)(數(shù)據(jù)來源于文獻(xiàn)[29])測試SCA反演實(shí)際數(shù)據(jù)的效果。原始地震記錄如圖7a所示。采集震源為錘擊震源,選用的檢波器主頻為4.5 Hz,共24道,道間距為1 m,偏移距為10 m。圖7b為提取的頻散能量,從圖上可看出高階波能量弱,基本不發(fā)育,反演選用基階波數(shù)據(jù)。反演時(shí)算法參數(shù)同前面理論模型測試時(shí)參數(shù)設(shè)置相同,模型參數(shù)按照文獻(xiàn)[30]進(jìn)行設(shè)置,所有地層的縱波速度及密度分別設(shè)置為1 440 m/s和1 850 kg/m3,其余參數(shù)如表5所示。表5中的層厚及S波速度數(shù)據(jù)來源于文獻(xiàn)[29]。
表5 Arnarb?lidi地區(qū)前人[29]反演模型參數(shù)以及搜索范圍[30]
圖7 Arnarb?lidi地區(qū)地震數(shù)據(jù)(a)與頻散能量圖(b)[29]
圖8a中可看出反演所得頻散曲線(黑點(diǎn))和實(shí)測頻散曲線(實(shí)線)擬合程度非常高;反演模型參數(shù)(圖8c中黑點(diǎn))在淺層與文獻(xiàn)[29]的參考模型(圖8c中實(shí)線)參數(shù)差別較小,只有第9層橫波速度出現(xiàn)較大差別,推測原因是基階頻散曲線低頻信息不夠充分,故對于深部結(jié)果約束性不強(qiáng)所導(dǎo)致。SCA迭代到60次時(shí)就基本收斂,其均方根誤差僅只有2.90 m/s,這小于文獻(xiàn)[29]中所得到的5.18 m/s,故將SCA應(yīng)用于實(shí)際瑞利波頻散曲線反演中的效果相較于前人結(jié)果具有更為優(yōu)異的表現(xiàn)。
a—反演所得頻散曲線;b—最小目標(biāo)函數(shù)值隨迭代次數(shù)變化情況;c—反演的橫波速度剖面
接著使用美國懷俄明地區(qū)數(shù)據(jù)(數(shù)據(jù)來源于文獻(xiàn)[31])進(jìn)一步測試SCA反演實(shí)際數(shù)據(jù)的效果。原始地震記錄如圖9a所示,采集地震數(shù)據(jù)時(shí),使用了48個(gè)檢波器,震源為錘擊震源,偏移距、道間距都為0.9 m。圖9b為提取的頻散能量,由于高階數(shù)據(jù)出現(xiàn)了模式接吻現(xiàn)象,文中僅選用基階波數(shù)據(jù)進(jìn)行反演。反演時(shí)算法參數(shù)同前面理論模型測試時(shí)參數(shù)設(shè)置相同,反演的模型參數(shù)如表6所示。
表6 懷俄明地區(qū)前人[32]反演模型參數(shù)及搜索范圍
圖9 懷俄明地區(qū)地震數(shù)據(jù)(a)與頻散能量(b)[32]
反演結(jié)果如圖10所示,從圖10a和10b中可看出,反演得到的頻散曲線(圖10a中的實(shí)線)與實(shí)測數(shù)據(jù)頻散曲線(圖10a中的點(diǎn)線)擬合得很好、精度較高,反演迭代20多次后就已經(jīng)收斂,收斂精度為4.40 m/s。在淺部,SCA重建的橫波速度模型(圖10c中的實(shí)線)與測井資料(圖10c中帶菱形的折線)都吻合得很好,在5 m往深后,重建模型依然能較好地貼近測井?dāng)?shù)據(jù)??芍猄CA用于實(shí)測數(shù)據(jù)反演具有較高可信度。
a—反演所得頻散曲線;b—最小目標(biāo)函數(shù)值隨迭代次數(shù)變化情況;c—測井?dāng)?shù)據(jù)與反演橫波速度剖面對比
本文將一種新的群智能優(yōu)化算法——SCA用于瑞利波頻散曲線反演的研究當(dāng)中。在進(jìn)行具體的反演工作時(shí),通過設(shè)置較大的模型參數(shù)搜索范圍以模擬無先驗(yàn)信息或者僅有少量已知信息的情況,并同時(shí)反演厚度、密度、橫波速度以及縱波速度以更好適應(yīng)真實(shí)情況(但只選取對頻散曲線有顯著影響的層厚和橫波速度來進(jìn)行分析與評價(jià))。利用多個(gè)理論模型(從簡單模型到復(fù)雜模型,從無噪聲數(shù)據(jù)到含噪聲數(shù)據(jù))對SCA性能從理論的角度進(jìn)行了詳盡的測試,最后再采用來自冰島Arnarb?lidi和美國懷俄明地區(qū)的實(shí)測數(shù)據(jù)檢驗(yàn)了SCA的實(shí)用性。測試結(jié)果表明:
1)基于正余弦函數(shù)組合,應(yīng)用多個(gè)隨機(jī)參數(shù)和自適應(yīng)變量調(diào)整尋優(yōu)過程的SCA對于頻散曲線反演搜索空間具有十分優(yōu)異的探索和開發(fā)能力,在獲得高精度地層模型參數(shù)的同時(shí)還保證了收斂速度以及較強(qiáng)的穩(wěn)定性,具有良好的發(fā)展?jié)摿Α?/p>
2)SCA相較于傳統(tǒng)的粒子群算法在瑞利波頻散曲線反演中的表現(xiàn)具有更高的精度、更快的收斂速度以及更強(qiáng)的穩(wěn)定性。