蚱蜢算法在瑞雷波頻散曲線反演中的應用

于東凱 宋先海*② 張學強 趙素濤 蔡 偉

(①中國地質大學(武漢)地球物理與空間信息學院，湖北武漢 430074；②中國地質大學(武漢)湖北省地球內部多尺度成像重點實驗室，湖北武漢 430074)

0 引言

瑞雷波由縱波與橫波的干涉形成，沿著自由表面傳播。英國學者瑞雷首先發現并證明均勻半空間中瑞雷面波的存在[1]。瑞雷波勘探具有無損、高效、經濟等特點，被視為未來最有希望的技術之一[2]，現已被廣泛用于淺地表、巖石圈橫波速度結構信息[3-5]及地層品質因子[6]的計算、路基路面檢測[7]、基巖面刻畫[8]等。與常規勘探方法相比，瑞雷波勘探具有分辨率高、應用范圍廣、檢測設備簡單、檢測速度快等優點[9]。近年來，諸多學者進行了相關研究，如瑞雷波反演[10-11]、瑞雷波頻散研究[12-14]、瑞雷波數值模擬研究[15-17]、地形對瑞雷波傳播的影響[18]等。

瑞雷波勘探可分為三個過程：面波數據采集、頻散曲線提取和頻散曲線反演[19-20]。拾取瑞雷波頻散曲線后，利用瑞雷波頻散特性可求得近地表橫波速度，是許多工程的關鍵參數之一[21-22]。多數局部線性化優化方法被用來反演頻散曲線，并在各種商業軟件中占主導地位，如最速下降法、L-M算法等。然而，與大多數其他地球物理問題一樣，瑞雷波頻散曲線反演是高度非線性、多參數、多極值的，局部線性化方法易收斂到局部極小值，反演結果在很大程度上取決于初始模型的選擇以及偏導數計算的準確性。全局最優化算法由于沒有上述限制，被越來越多地應用于頻散曲線反演，如Yamanaka等[23]借助遺傳算法反演頻散曲線研究地下結構信息； Beaty等[24]利用模擬退火算法反演多模式瑞雷波頻散曲線； Shirazi等[25]通過人工神經網絡算法反演頻散曲線以獲得路面結構； Song等[26-27]先后將模式識別算法、粒子群算法應用于頻散曲線反演。傳統算法，如遺傳算法，在解決反演問題中具有較高的應用前景，然而往往存在計算量偏大、效率不高并且容易早熟收斂等缺陷，在反演中受到限制。

隨著智能算法的發展，新型算法不斷被應用于地球物理反演[28-31]。本文針對傳統算法中出現的問題，引入一種新型智能算法——蚱蜢算法，對瑞雷波頻散曲線進行反演以獲取近地表橫波速度剖面。對不同地質模型的理論頻散曲線進行反演試驗，并進行多次、多角度、深層次的探索，最后利用實測數據檢驗其實用性。

1 蚱蜢算法基本原理及驗證

1.1 基本原理

受蚱蜢不同時期覓食行為的啟發，澳洲學者Saremi等[32]提出了一種新的仿生學算法——蚱蜢算法(Grasshopper Optimization Algorithm，GOA)，算法模擬了蚱蜢在幼蟲和成蟲兩個階段搜索食物時所表現出的不同行為。蚱蜢是自然界中一種常見的昆蟲，它們成群結隊地出現，在幼蟲時依靠強有力的雙腿進行跳躍和移動，雖移動緩慢，但能夠較為詳盡地對所經區域是否存在食物進行判斷；當其進入成蟲期后長出翅膀，可在空中飛翔，移動變得十分迅速，且移動范圍迅速擴大，搜索范圍得到大幅度的擴大，但同時對食物的搜索精確度會得到一定程度的削弱。

蚱蜢算法將蚱蜢個體模擬為搜索算子，將成年蚱蜢大范圍且迅速移動和幼蟲期局部移動模擬為算法尋優過程中隨迭代次數的增加所進行的全局探查與局部搜索兩個階段。在全局探查階段，算子在模型內進行大范圍的移動，獲取全局信息，并最終固定為某局部區域；在局部搜索階段，算子的搜索范圍縮減至某局部區域，在該區域內不斷迭代優化所得解，不斷提高解的精度。算法主要包括兩個要素：自身所處位置以及與其他算子的位置關系。其中，自身位置體現了算子當前的優劣，蚱蜢與其他算子的位置關系可表示為：搜索算子的移動受群體中其他算子的影響，當算子相隔較近(距離小于2.079個單位)時會出現彼此排斥現象，此時稱蚱蜢處于排斥域范圍內；當算子相隔大于該范圍時又會出現吸引現象，稱其位于吸引域內；當距離恰好相差2.079個單位時，彼此之間不會發生吸引或排斥現象，蚱蜢位于舒適域。算法通過不斷調整自身位置及與其他算子的相對位置實現優化。需要說明的是，這兩個階段之間沒有明顯的界限，隨著搜素區間逐步縮小，算法由全局階段逐漸轉至局部階段。蚱蜢算法的優化原理描述如下。

蚱蜢算子xi受蚱蜢算子xj的影響因子定義為

(1)

假設i個蚱蜢算子當前位置為xi,k， 移動時會受到所有其他算子的影響，因此移動后的位置為

]}+xi,k

(2)

式中：i和j表示不同蚱蜢算子；NL表示算子總數；di,j表示第i個與第j個算子的空間距離；h和g均為評價算子受其他算子影響的參數，前者表示吸引力強度，后者是衡量空間距離的尺度，二者與影響因子F(di,j)緊密相關，不同組合會導致算子不同的行為習慣，算法中g通常取1.5，h取0.5；u(k)和l(k)分別表示算子的第k個分量在區間內的上、下限；xi,k表示當前蚱蜢算子所處位置；c1和c2為自適應的遞減變量，c1與粒子群算法中慣性權重相似，作用是減少蚱蜢算子圍繞目標值的移動，c2作用是收縮算子群體的舒適域、排斥域和吸引域，通過對自身的調節可在有限迭代次數中達到全局探查以及局部搜索的相對平衡。

總體來說，式(2)的前半部分實現了其他蚱蜢算子的約束作用，后半部分則對當前算子移動至全局最優解的行為趨勢起到一定的控制作用。c1和c2雖在算法中所起作用不盡相同，但均為隨迭代進行自適應衰減的變量，可統一表示如下

(3)

式中：Iiter代表當前迭代次數；NI為最大迭代次數；cmax為自適應參數最大值，本文cmax取值為1；cmin為最小值，文中取1.0×10-5。可根據解決問題的不同設定變量c1或c2的起始值。

本算法主要受三個參數控制：種群數量NL、迭代次數NI以及自適應參數c1和c2。種群數量與迭代次數對算法的影響相似：隨著種群數量取值的增大，算法在變量空間中搜索會變得詳盡，但同時也可能出現過多的冗余搜索，影響算法的效率；若種群數量取值過小，則會出現搜索不充分，無法找到全局最優解。因此，取值要根據所解決的具體問題來設定。自適應參數的取值與算法的尋優功能緊密相關，具體來說，為避免所得解陷入局部最優，設置如下機制：蚱蜢群體對某算子會產生排斥或吸引現象，其中對算子的排斥作用能有效避免局部最優解的產生，阻止算法陷入局部最優以及停滯不前。同時，由于自適應參數的引入，算法會縮小與迭代次數成比例的排斥面積，并不會對算法的收斂性造成過大影響，因此蚱蜢算子在開始階段能有效避免局部最優解，在最后的局部搜索階段也能改進優化所得解。

蚱蜢優化算法的基本流程如下。

(1)參數初始化。在模型搜索空間內隨機生成蚱蜢算子的初始位置xi(i=1,2,…,NL)，每個算子的位置可表示為xi=(xi,1,xi，2，…,xi，k,…,xi，D)，這里D表示算子中包含的變量個數。

(2)根據蚱蜢算子所處位置，計算各個算子的目標函數值，并找出最優解及相應的算子位置。

(3)更新變量c1和c2，計算蚱蜢之間空間距離，并對所有蚱蜢算子遍歷，根據式(2)更新空間位置； 同時，若新位置躍出上下邊界，則給予約束。

(4)當達到最大搜索次數或解滿足給定條件時，轉至步驟(5)； 否則，搜索次數增加1，轉至步驟(2)，進行下一次搜索。

(5)輸出全局最優解及其對應的最優擬合值。

1.2 算法測試

不同算法在處理問題上的能力存在差異，同一算法處理不同問題的能力也有高低。故在全局優化算法領域，通常使用一些特定的函數對算法進行測試，基于生成解的精度對算法的優化能力進行評價。本文主要研究將蚱蜢算法應用于多極值的頻散曲線反演，多模態復合測試函數與反演中的目標函數存在一定的相似性，即具有多個局部最優解，因此選取兩個含有多峰值的函數進行全局尋優能力測試。

(1)Griewank函數。如圖1a所示(為方便顯示，僅顯示局部結果)，該函數存在許多局部極小值，其數目與問題的維數有關，全局最小值0在(x1,x2,…,xD)=(0,0,…,0)處，因此該函數是典型的非線性多模態函數，具有廣泛的搜索空間，通常被認為是優化算法很難處理的復雜多模態問題，在程序中搜索算子個數NL取500。其函數表達式為

xi∈[-600,600]

(4)

(2)Rastrigrin函數。如圖1b所示，該函數是一個多峰值函數，在(x1,x2,…,xD)=(0,0,…,0)處取得全局最小值0， 在{xi∈(-5.12, 5.12),i=1,2,…,D}范圍內大約有10D個局部極小值，此函數與Griewank函數類似，也是一種典型的非線性多模態函數，峰形呈高低起伏不定跳躍性地出現，所以很難優化查找到全局最優解。其函數表達式為

xi∈[-5.12, 5.12]

(5)

測試中，種群數量NL取100，迭代次數NI取700。表1列出了蚱蜢算法對兩個測試函數的優化結果，該結果為5次試驗的平均值，可看出算法的所得解均收斂到函數最小值附近，收斂誤差接近于零，說明算法具備較強的全局尋優能力。圖2a為Griewank函數的算法擬合誤差隨迭代次數的變化，可見迭代逐漸收斂并最終趨近于零，說明GOA法通過多次迭代能夠有效提高解的精度。另外，圖中也顯示出擬合誤差發生了不確定性變化，一方面是因為算法中引入了貪婪準則，保留當前最優解，算法在迭代之后所得解的精度未得到提高時會出現誤差不變的情況； 另一方面，高排斥率導致蚱蜢算子在全局范圍內發生大范圍的移動，使迭代誤差發生較大的改變。圖2b所示為算法在多峰值、多模態測試函數Rastrigrin函數中的表現，圖中每個點均代表一個蚱蜢算子在某次迭代中的搜索位置。可以看出：搜索過程中算子已基本覆蓋所有目標區域，可見算法具有較高的搜索能力；圖中紅點表示全局最優解，算子在該點附近大量聚集，說明算法并未陷入局部最優解，而是逐步逼近最優解。另外，圖中有“抱團”現象，即算子在某些區域聚集，原因可能是該點為局部最優點，蚱蜢算子在該處附近迭代次數較其他位置稍大。

圖1 測試函數圖

圖2 蚱蜢算法中Griewank函數的迭代誤差曲線(a)和Rastrigrin函數的算子位置圖(b)

測試函數Griewank函數Rastrigrin函數理論極值00實際極值0.07880.0037算法所得解(4.17×10-4,7.01×10-3)(2.30×10-6,-1.35×10-6)

2 蚱蜢算法在面波反演中的實現

在瑞雷波勘探的應用中，依據面波頻散資料推算工程場地的剪切波速度剖面時，通常認為地下為層狀介質。同時，基階波在面波頻散資料中能量最高、最易觀測，也應用最為廣泛，故本文主要采用基階模式頻散曲線反演獲取地下層狀地質結構。夏江海等[6]的研究表明，與瑞雷波頻散曲線特征關系最密切的是地層的橫波速度(vS)和厚度(H)，這兩個變量可通過反演求取，地層密度、泊松比為已知參數(對頻散曲線影響微弱，可通過工區信息大致估算)。另外，為更符合實際情況(通常沒有足夠的先驗信息估計淺地表的橫波速度和厚度準確范圍)，文中使用較大的搜索范圍，在下文理論模型測試中，搜索區域的下限和上限設定為與真實值相差50%或更多。如前文所述，蚱蜢算法主要受到三個參數控制：NL、NI和c。反演程序中NL取值為10D；NI取值200；cmax和cmin分別取1和1.0×10-5。

頻散曲線反演是一個求解目標函數最小值的優化過程，算法的好壞取決于能否在給定范圍內充分搜索以找到反演參數的最優解。同時，由于反演中存在局部極值的干擾，算法能否跳出局部極值、找到全局最優解也是評價其質量的標準之一。目標函數是通過算法反演得到的模型能否精確解釋觀測資料所確定的，應用最為廣泛的均方差函數是

(6)

(7)

目標函數具有以下幾個特點：

①目標函數為M個頻點(即頻散曲線上的采樣點)上的相速度共同作用。值得注意的是，采樣數目過大并不會對反演過程產生顯著影響，反而會影響反演效率[26]。由于所研究面波的頻散曲線頻率多在5～100Hz范圍，試驗發現頻點數取25～30能平衡算法的高效性與采樣的保真度。

針對目標函數的特點，蚱蜢算法在處理反演問題上具有獨特優勢，體現在算法能夠保證目標函數不困于局部極值，而是進行全局范圍的搜索；同時不斷優化所得解，提高反演精度。算法引入排斥域、吸引域和舒適域機制：當某位置上多個算子之間距離過小時(在反演中表現為迭代誤差相近)會相互排斥，算子在下次迭代中遠離該位置，在算法初期需深入分析變量空間，避免受困于局部極值；當距離過大時認定進入吸引域，算子之間會產生吸引力，體現在中后期的“開發工作”中，隨著算法不斷運行，大量的非可行區間被排除，搜索范圍會不斷被壓縮至最優解附近，算法對此區域不斷開發，對所得解不斷優化以求提高反演精度；當距離適當時，算子彼此之間不會吸引或排斥，此時位于舒適域。為更好地平衡“開發”(即全局探查)與“改進”(即局部搜索)兩個過程，本算法配備了一個適應性縮減三個作用域面積的系數，隨著迭代的進行，算法會逐漸由全局探查轉為針對最優解附近區域的局部搜索。由此可知，通過群體獲得的最佳解決方案是蚱蜢開發和改進的結果。因此，在解決一系列反演問題時，蚱蜢算法具有超越當前幾種算法的潛力。

2.1 理論模型反演

為檢驗蚱蜢算法在瑞雷波頻散曲線反演中效果，使用兩個理論模型進行測試，二者均為實際淺層工程勘察中經常遇到的地質模型。如表2所示，模型A代表橫波速度遞增的四層地質模型；模型B代表在兩個較硬地層之間夾雜低速軟弱地層，該地層剪切波速度低于上覆地層速度，常見的有路面結構等。

對理論模型進行蚱蜢算法反演試驗，模型A的數值模擬及反演結果見圖3和圖4，模型B的數值模擬及反演結果見5和圖6。圖3a、圖5a分別為兩個模型波場正演得到的地震記錄，以此模擬實測資料并作為頻散曲線反演的依據。圖4a、圖6a為擬合頻散曲線，圖3b、圖5b是通過高分辨率線性拉東變換提取的模擬地震記錄得到的頻散能量圖。值得注意的是，由于含低速軟弱夾層模型中頻散曲線的特殊性，即基階模式在高頻范圍內(該模型中為頻率大于60Hz時)攜帶的能量并不是最高的，因此文中拾取的基階波頻散曲線也有一定的范圍。

表2 兩種地質模型參數及反演搜索范圍

圖3 模型A模擬地震記錄(a)及頻散能量圖(b)

圖4 基于蚱蜢算法的模型A反演結果

圖5 模型B模擬地震記錄(a)以及頻散能量圖(b)

圖6 基于蚱蜢算法的模型B反演結果

可以看出，在沒有足夠先驗信息的情況下，頻散曲線吻合度較高，且模型中的參數均可被精確地反演和重建。就解的精確性而言，兩個模型中參數的最大相對誤差相對分別為2.50%、2.89%。由此可見，通過蚱蜢算法獲得的地下橫波速度與真實情況相差不大，說明本算法可以對理論瑞雷波頻散曲線進行很好地反演，對頻散曲線影響較大的參數(各層橫波速度以及層厚度)均能被精確地反演和重建。

2.1.1 頻散曲線聯合反演測試

瑞雷波是頻散且多階的，某頻率下的最小傳播速度被稱為基階模式，其余按相速度大小依次稱為一階高階模式、二階高階模式……。在瑞雷波頻散曲線的應用中，由于基階波能量較強且易觀測等原因，通常使用基階波頻散曲線進行非線性反演。但在某些地形條件下，如文中的模型B，高階模式的頻散曲線在高頻范圍內攜帶了更多的能量，此時若只使用基階波進行反演可能存在有效信息不足的問題。為進一步驗證算法的反演能力，本次測試針對模型B，利用蚱蜢算法對基階和高階模式頻散曲線進行聯合反演。聯合反演的頻散曲線共N條，故目標函數相應變為

(8)

圖7b所示為橫波速度剖面的對比。表3為基階波反演與聯合反演的對比，可以看出聯合反演結果較基階波反演得到進一步優化，各參數基本接近真實值，反演誤差趨于零，此結果進一步證實蚱蜢算法反演精度較高。

圖7 基于蚱蜢算法的模型B聯合反演結果

反演方式vS1m·s-1vS2m·s-1vS3m·s-1vS4m·s-1vS5m·s-1H1mH2mH3mH4m基階波反演195.48147.14204.30304.01398.572.012.013.996.01聯合反演 199.27149.58201.29301.74399.502.002.014.006.01

2.1.2 含噪測試

由于野外環境復雜，采集到的地震記錄中不可避免含有噪聲，含噪數據可能降低算法精度及穩定性，并影響反演結果。對一種新的反演算法需檢驗其抗噪能力。本文對理論頻散曲線進行加噪處理(加入5%、10%及20%噪聲)，即對各頻率下的相速度以理論值為基準在一定區間內進行隨機程度的增減，以模擬實際情況中混雜了環境噪聲的相速度，含噪頻散曲線處理結果見圖8。

圖8a中雜亂不規則的藍色、藍綠色以及黑色實線所示為加入不同程度隨機噪聲的頻散曲線，圖8b為算法對含噪數據的反演結果。由圖8可見，橫波速度與真實模型之間有一定的偏離，三種不同噪聲水平的數據反演結果最大相對誤差分別增加到3.5%、9.5%和13.0%，由此可見噪聲對反演解的準確性有一定的影響，但反演速度剖面均未實質性地偏離真實模型，因此由算法得到的反演結果仍然具有較高的可信度，說明蚱蜢算法對隨機噪聲具有一定的容忍度，用于反演瑞雷波頻散曲線具有較高的穩定性。

2.1.3 搜索區間測試

前文設定的搜索范圍是根據前期測井等數據對整個工區層位進行推斷，并以此估算地震記錄下對應各層橫波速度及厚度的大致范圍。倘若地質環境復雜，地層可能發生突變，則上述層位信息在當前地震記錄下是不適用的，有必要設定更多數量的薄層以求準確模擬近地表地質情形。這是因為應用快速矢量傳遞算法計算頻散曲線時需要輸入各層橫波速度和厚度，比如n層模型需輸入(2n-1)個參數，故在反演之前有必要對模型層數進行設定。本文為進一步檢測算法的反演能力，以模型A為例，將模型層數設定為7，各層層厚均設定為2m，且橫波速度搜索范圍均設定為100～800m/s，頻散曲線反演是在給定的范圍內找出最優模型，因此要預先設定反演中橫波速度的搜索范圍，反演結果如圖9所示。由圖可見，在先驗信息有限、搜索空間更為廣泛的前提下，反演過程耗時增加，反演結果仍十分接近真實值，并未發生顯著變化，說明算法的全局搜索能力較強，反演結果具備較高的可信度。

圖8 含不同噪聲水平的模型A理論數據反演結果

圖9 擴大搜索空間時基于蚱蜢算法的模型A反演結果

2.1.4 頻散曲線缺失頻段測試

在實地勘測中獲得的瑞雷波頻散曲線具有先天的不完整性，存在部分頻段受干擾或缺失等情況。實際上，在松散的沉積物上施加低頻源(如大錘)所產生的頻譜往往會缺乏高頻率數據。針對丟失頻段的瑞雷波頻散曲線反演時，首先應考慮目標層的敏感頻率范圍。本文針對模型A首先分析基階頻散曲線對各層橫波速度的敏感度。從圖10a可以看出，基階波對各層橫波速度變化敏感的頻段隨深度變化而變化，表層的敏感頻段為高頻，由于穿透深度的影響，深部層敏感頻段的頻率變得更低更窄，即面波的高頻成分對淺部地層的橫波速度較敏感，而低頻部分對深度層和半空間的橫波速度更敏感。基于此，截取0～50Hz范圍的頻散曲線進行基階波頻散曲線反演，通過對比橫波速度剖面分析蚱蜢算法反演能力。從圖10b可見，在缺失部分頻段的情況下，多次試驗結果雖然出現一定程度的波動，誤差較之前也增大，但均未實質性地偏離真實值，且多次反演的均值與真實剖面吻合度較高，說明算法的反演能力值得信賴。

2.2 反演結果分析

通過理論模型運算及四種測試全方位地分析了蚱蜢算法的反演能力，其中穩定性、收斂性和有效性是評價反演算法的重要指標。為避免單次反演的偶然性，對上述測試均進行20次獨立反演運算，并取均值討論。受篇幅限制，只列出速度遞增模型A的反演結果統計表(表4)。可以看出：對于無噪數據，模型反演結果與真實值相符，相對誤差和均方根誤差均趨近于零，說明結果準確可靠； 對于含噪測試中含噪數據(10%隨機噪聲)，模型的7個參數中相對誤差最大不超過9.5%，均方差最大不超過14.32，這說明噪聲會對反演精度造成一定的影響，但所得解與真實解相差并不大，精度能滿足實際需求。原因在于蚱蜢算法引入舒適域、排斥域和吸引域三個作用域，算子會受到群體中其余所有算子的共同作用，全局尋優能力更強，不易出現早熟收斂，能夠有效避免其他算法中由于隨機化個體單一、受當前最優解影響過大導致的全局尋優能力不足等問題。

圖10 基于蚱蜢算法的模型A反演結果(缺失頻段)

參數真實值不含噪聲含10%隨機噪聲反演均值相對誤差/%標準差反演均值相對誤差/%標準差vS1/(m·s-1)200.00203.841.921.57208.114.0614.32vS2/(m·s-1)250.00252.400.963.02257.402.965.60vS3/(m·s-1)300.00303.241.082.14290.243.2510.50vS4/(m·s-1)400.00396.810.803.26394.011.507.86H1/(m)2.001.952.501.202.042.001.02H2/(m)2.002.031.500.372.199.500.79H3/(m)4.004.061.500.414.082.000.85

另外，由搜索區間測試中的收斂曲線可以看出，無論搜索區間的上、下邊界是否靠近真實值，算法的最優擬合值在初期的迭代過程中均迅速下降，說明算法在此階段迅速排除搜索空間內大量非可行區間，類比于蚱蜢在成蟲期的覓食中能夠迅速排除非食物源區域。在隨后迭代過程中誤差下降趨勢變緩，且在后期趨近于零，說明算法的反演過程在此之前已基本完成，各個參數已被精準地推算。值得注意的是，任何一種群智能算法都會出現多次迭代中解的精度沒有得到提高的情況，這是因為在某幾次或多次迭代過程中，種群中沒有找到更優個體，根據貪婪準則保留當前最優解，解的精度沒有得到提高，故而誤差不變。

同時，頻散曲線反演是一個非線性、多極值問題，不同的模型參數組合也可能獲得相似的擬合值，因此還需要分析每次反演所得解中模型參數的分布概率，以驗證反演結果的有效性。圖11所示為模型A中各個參數在真實值附近的分布概率。圖11a所示為第一層橫波速度vS1在20次反演結果中的分布概率，真實值為200m/s，反演所得解在真實值±10m/s范圍內的概率接近0.80。由圖11可見，其余模型參數在其各自真實值附近(速度±10m/s，厚度±0.1m)的分布概率分別為0.75、0.60、0.70、0.65、0.70和0.65。可見反演所得模型參數均大致分布在真實值附近，進一步說明蚱蜢算法是一種可靠、準確、收斂、有效的反演算法，能夠有效地應用于瑞雷波頻散曲線的反演。

圖11 模型A反演參數的分布概率直方圖

3 與粒子群算法的對比

為驗證蚱蜢算法在處理地球物理反演問題中的適用性，將該算法與粒子群算法(Particle Swarm Optimization，PSO)進行比較。粒子群算法是由Kennedy等[34]提出，后得到眾多學者的改進和優化[35]，已成為當今最流行、最成熟的自然啟發優化算法之一，并在地震勘探中已獲得廣泛應用[36-38]。前人已成功將其應用于頻散曲線反演[27]，本文就基于新型算法與粒子群算法的頻散曲線反演結果進行對比，證明前者具有較高的可信度。

這兩種算法具備如下幾個相同點：均為全局優化算法，算法在全局解空間內進行搜索，并將搜索重點集中在性能高的部分；都是隨機搜索算法，且具備記憶性能，最優解都被保留下來；具有隱含并行性搜索特性，搜索過程均從解空間內的一個集合開始。但蚱蜢算法較粒子群算法而言，在處理地球物理反演問題上具有不易陷入局部極值的優勢。后者根據粒子本身歷史最優位置和當前群體內的最優算子的位置決定下一步的移動方向，因此可能會出現如下現象：由于目標函數存在大量的局部極小值，若最優算子進入某個極值位置，則算法接下來會在該位置附近多次移動，反演效率不高。而蚱蜢算法中算子會受到群體內其他所有算子的共同作用，體現在程序中算子每一步的移動方向不僅受當前最優算子的影響，也會受到群體內其余算子的共同作用，一定程度上降低了上述現象產生的可能性。同時，當算法運行至中后期時，算法通過縮小與迭代次數成比例的排斥面積，會減少排斥現象的產生，故不會對算法的收斂性造成過大的影響。

以模型B為例，分別采用上述兩種算法對模型B不含噪聲的理論數據進行10次獨立反演測試。為排除其他因素的干擾，采用相同的群體規模、搜索空間和迭代次數，M取值30。測試分析一次反演的結果而不是10次反演迭代誤差的均值，目的是更準確直觀地分析算法的性能。反演結果如圖12所示，由圖可以看出，如果迭代截至50～100次(圖12a)，迭代誤差的整體變化情形與算法優化所得反演解的能力基本一致，兩種算法效果不相上下，說明這兩種算法均具備較強的收斂能力和收縮搜索區間性能；但迭代進入中后期(圖12b)，即兩種算法運行100～200次，蚱蜢算法比粒子群算法收斂速度更快，且收斂精度也有一定提高，說明蚱蜢算法在處理地球物理反演問題時具備較高的尋優性能。

圖12 模型B不含噪聲時的GOA與PSO迭代誤差對比

4 實測資料反演

前文充分驗證了蚱蜢算法在理論瑞雷波頻散曲線反演中可行且有效，所得解精度高。為進一步檢驗算法的適用性，利用蚱蜢算法對美國懷俄明某地獲得的地震數據進行了分析。圖13a所示為采集的地震記錄，數據采集采用了48個8Hz垂直分量檢波器，道間距為0.9m，最小炮檢距為0.9m，震源采用錘擊震源。圖13b所示為對應的頻散能量圖，通過采集能量最大點提取基階波頻散曲線。與反演理論模型類似，固定地層的密度及泊松比，只反演橫波速度vS與厚度H兩個參量。根據測井數據將地層大致劃分為5個層位，并粗略估計搜索范圍。反演模型的搜索空間、泊松比以及密度等參數設置見表5。

圖14a所示為算法所得解與實測數據的基階波頻散曲線擬合情況，通過拾取頻散能量圖中能量最大值(圖13b中圓點)得到不同頻率下的相速度，可以看出頻散曲線擬合情況較好。文中采用基階波反演的原因在于能量圖中出現較為嚴重的“模式接吻”以及“模式跳躍”現象，在拾取高階頻散曲線的過程中很容易出現模式誤判，從而影響反演結果，若采用基階波反演精度會更高。圖14b為兩種算法反演時迭代擬合誤差隨迭代次數的變化曲線，為保證反演的精準性、降低因單次反演中出現的偶然誤差對結果造成影響，針對野外數據進行了多次反演，舍棄異常的反演解，并對剩余解進行均值化處理，因此迭代誤差曲線較圖12更為平滑。圖14c為橫波速度剖面對比，可見通過有限頻段的頻散曲線反演所得速度模型與測井資料吻合較好，尤其是淺地表部分。相較于PSO算法，GOA算法能夠準確地分辨厚度為3～5m的軟弱夾層。值得注意的是，較前文測試，此次實測資料反演的信息量有限，與模型試算結果相比反演效果并不突出，但借助有限的頻散信息所得的反演結果能基本滿足實際工作需求。另外，較深地層反演速度剖面與測井資料出現了一定偏離，這是由于采用錘擊震源，瑞雷波頻散能量在低頻段分辨率低，導致提取誤差增大；另外一個原因是算法只設定5個層位，較深的地層被認定為均勻介質，即橫波速度沒有變化。

表5 GOA反演中模型搜索范圍及模型參數

圖13 懷俄明某地區地震記錄(a)與頻散能量圖(b)

圖14 懷俄明某地區地震數據蚱蜢算法反演結果

5 結論

本文介紹了一種全新的仿生學算法——蚱蜢算法，并通過兩個測試函數驗證其全局搜索能力，證明了其優異的收斂能力以及全局尋優能力。將蚱蜢算法應用于瑞雷波頻散曲線反演，通過設置較大的搜索范圍以模擬實際工作中缺少足夠先驗信息的情況，并進行多角度、深層次的探索試驗，深入分析了算法的反演能力。最后通過實測數據檢測了本算法的適用性。理論數據和實測資料反演的結果表明：

(1)蚱蜢算法具有較高的全局尋優能力以及抗局部最優解能力，所得解與函數最小解十分接近，說明這是一種有效的全局優化算法；

(2)蚱蜢算法可應用于瑞雷波頻散曲線反演，與傳統算法相比，反演過程中目標函數曲線收斂更加穩定，最終迭代誤差更接近于零。證明基于蚱蜢算法的瑞雷波頻散曲線反演具有較高的精度，是一種簡單、穩健的反演算法。

同時，蚱蜢算法也存在一定的缺陷，比如算法的收斂過程相對較長。在今后的研究中，可與其他算法如遺傳算法、粒子群算法或人工神經網絡算法相結合，探討如何進一步提高算法的收斂效率。